# FindRoot inside Table

(MWE at the end of the post)

I need to solve a non-linear equation $f(y;x_1,x_2,..,x_5)$ in one variable $y$ and then compute 4 new output expressions, for over 60 different initial parameter inputs of $x_i$. The 4 output variables which are as follows:

$g_m(y,x_1,..,x_5) \; \forall m$

The main equation will solve for the variable $y(i)$:

$f(y(i),x_1(i),x_2(i),...,x_5(i))=0 \; \forall \; i \in \{1,60\}$

Now my

$f(y(i),x_1(i),x_2(i),...,x_5(i))=A1(y,x_1,..,x_5)*y^{k_1} + B1(y,x_1,..,x_5)*y^{k_2}-c1(x_1,..,x_5)$

I define variables

A1=A1(y,x1,..,x5);
B1=B1(y,x1,..,x5);
c1=c1(x1,..,x5);
g1=g1(y,x1,..,x5);
g2=g2(y,x1,..,x5);
g3=g3(y,x1,..,x5);
g4=g4(x1,..,x5);


I have data for $x_1(i),..,x_5(i)$ in a CSV file that I Import and Table straight into the variable names making them into lists.

datatemp = Import["C:\\Documents\\2012U26G0.csv"];
j = Dimensions[datatemp][[1]]
kk=2
x1 = Table[datatemp[[i, 3]], {i, kk, j}]
x2 = Table[datatemp[[i, 4]], {i, kk, j}]
x3 = Table[datatemp[[i, 5]], {i, kk, j}]
x4 = Table[datatemp[[i, 6]], {i, kk, j}]
x5 = Table[datatemp[[i, 7]], {i, kk, j}]


I think this automatically makes my earlier defined formulae for $A1$,$B1$ and $g_m$ into a list of formulae with the only unknown being $y$ and it makes $c1$ into a list of constants since $c1$ was only dependent on $x_i$.

Now, what I would like to be able to do is the following. Give some initial search point for my FindRoot for $i=1$

sol={1,1,1,1,120}


As you will see in a second, I only care about

sol[[5]]


Due to continuity, the roots move monotonically with $i$ so once I find one root, I can get a sense of where to look for the next one so I substitute the previous solution into the search for the next one. Also, in one shot I compute the 4 output variables I need. So when the Table is run, in one shot I have all the output data I want.

outputdata=Table[sol={g1,g2,g3,g4,y} /. FindRoot[ A1[[i]]*(y[[i]])^(k1) + B1[[i]]*(y[[i]])^(k2)==c[[i]] , {y, sol[[5]]+10, sol[[5]], sol[[5]]+20}], {i,1, 60}]


This process worked for a charm for little while but for a certain parameter space (by parameter I don't mean the $x_i$ I used earlier but a host of $\gamma$s and$\beta$s in my equations that I have suppressed so far), it has started giving me errors, 1/0 infinity type stuff, because of some assignment issues. Is there a clean and correct/good way to do this? I wanna be able to import a ton of data, Table my findroot to compute a whole bunch of data and Export it real fast. Please please please help!

MWE

f = y^(3.1276)*(A1) + y^(-0.5875)*(B1) + (c1)^2;
A1 = x1/y + x2*y + 3*x3;
B1 = x1*x3 + 1/(y*x2);
c1 = x1^3 + x2^5 - x3;
g1 = y^(x1) - x3*x2;
g2 = x1/y;


Imported Data Below:

x1 = {89, 88, 87}
x2 = {0.048334203, 0.048515211, 0.048707816}
x3 = {-19486.2273, -19742.04035, -20016.22863}


When I do this, I can see what the curves look like:

 Plot[Table[y^(3.1276)*(A1[[i]]) + y^(-0.5875)*(B1[[i]]) + (c1[[i]])^2, {i, 1,
3}], {y, 150, 180}]


This is what I want to do. To be able to Table a whole bunch of output in one shot:

dataoutput=Table[{g1[[i]], g2[[i]], y} /. FindRoot[y^(3.1276)*(A1[[i]]) +
y^(-0.5875)*(B1[[i]]) + (c1[[i]])^2, {y, 165, 150, 170}], {i, 1,
3}]


This is the result that I get, with some errors about accuracygoal and precisiongoal.

{{1623.03, -6.88842*10^19, 167.181}, {1530.37, -5.38632*10^19,
163.049}, {1431., -4.18019*10^19, 158.952}}


My MWE is working, just like my actual problem worked for a certain parameter space but now I am running into trouble. Is there a good way to Import a ton of Data, Table my FindRoot and generate a ton of output using Table and then Export my results? Thanks.

-
Please including a minimal example in Mathematica code. This is too vague to follow and all the "action" is inside the mystery file that is being imported. –  rm -rf Sep 9 '12 at 20:25
Thanks RM. I will post a MWE. There is nothing mysterious in the imported file but data. Just 60 values each for $x_1,x_2,...,x_5$ –  Amatya Sep 9 '12 at 20:37
You can reduce the complexity till you can settle on something small enough that demonstrates your problem. Like 3-4 values instead of 60, etc. (unless if the number 60 is the issue)... –  rm -rf Sep 9 '12 at 20:45
MWE included with 3 values of Imported data. –  Amatya Sep 9 '12 at 21:19

## Error Free:

You have almost solved your problem and I will try to help you streamline the whole process involving "Tons of data"! Lets first figure out the RHS of your equation in a symbolic form.

f = y^(3.1276)*(A1) + y^(-0.5875)*(B1) + (c1)^2;
A1 = x1/y + x2*y + 3*x3;r
B1 = x1*x3 + 1/(y*x2);
c1 = x1^3 + x2^5 - x3;
g1 = y^(x1) - x3*x2;
g2 = x1/y;
f // FullSimplify


(x1^3 + x2^5 - x3)^2 + (1 + x1*x2*x3*y)/(x2*y^1.5875) + y^2.1276*(x1 + y*(3*x3 + x2*y))

Now we set the parameter values.

x1 = {89, 88, 87};
x2 = {0.048334203, 0.048515211, 0.048707816};
x3 = {-19486.2273, -19742.04035, -20016.22863};
rhs = Table[y^(3.1276)*(A1[[i]]) + y^(-0.5875)*(B1[[i]]) + (c1[[i]])^2, {i, 1,3}];
scale = rhs /. y -> 150


{1.50954*10^11, 1.12911*10^11, 7.63367*10^10}

Now the trick to make FindRoot free of error message. We know for a an equation $f(x)=0$ one can scale the RHS by a constant say $\alpha \in \mathbb{R}$ to get $\frac{1}{\alpha}f(x)=\tilde{f}(x)=0$.We will scale each of your three equations by the above scaling factor respectively and then call FindRoot.

Plot[Evaluate[rhs/scale], {y, 150, 180}, Frame -> True,AxesStyle -> Dashed ]


MMA does not spill any error message!

eqs = (# == 0) & /@ (rhs/scale);
res = FindRoot[#, {y, 150}] & /@ eqs


{{y -> 167.181}, {y -> 163.049}, {y -> 158.952}}

However I was not able to reproduce your end result involving FindRoot. Notice that first coordinates are extremely big numbers due to the definition of g1.

dataoutput = Table[{g1[[i]], g2[[i]], y} /. res[[i]], {i, 1, 3}]


{{7.30459*10^197, 0.532358, 167.181}, {4.83167*10^194, 0.539714, 163.049}, {3.23723*10^191, 0.547335, 158.952}}

## Export and large data set:

To mimic your real data here I form a larger sample data of size $10000$. I assume that each parameter x1,x2,x3 follows a NormalDistribution with certain mean (here mean of x1,x2,x3 you supplied) and some arbitrary variance. Then I draw random sample from those distributions.

dataX1 = Floor@RandomVariate[NormalDistribution[Mean[x1], 3], 10^4];
dataX2 = RandomVariate[NormalDistribution[Mean[x2], .5], 10^4];
dataX3 = RandomVariate[NormalDistribution[Mean[x3], 16], 10^4];
TonsOfData = Transpose@{dataX1, dataX2, dataX3};


Now we form a function that does all the calculations and Apply it on TonsOfData. It took just arounf 7.5 seconds in my laptop.

With[{x1 = #[[1]], x2 = #[[2]], x3 = #[[3]]},
Rhs = (x1^3 + x2^5 - x3)^2 + 1 + x1 x2 x3 y/(x2 y^1.5875) +
y^2.1276 (x1 + y (3 x3 + x2 y));
Scaler = Rhs /. y -> 150;
sol = FindRoot[Rhs/Scaler, {y, 150}];
{y^(x1) - x3*x2, x1/y, y} /. sol
] & /@ TonsOfData; // AbsoluteTiming


{7.4364254, Null}

See how the roots are found gradually by FindRoot using Monitor for a data set of $1000$.

Now we simply export the result into a file called data.csv. Read the documentation to check which other data formats are supported by MMA.

Export["C:\\Users\\MMA\\Desktop\\data.csv", TonsOfResult, "Data"];


BR

-
PlatoManiac, can I have your babies? –  Amatya Sep 10 '12 at 19:53
Also, if it's not too much trouble, can you please write in the comments the code that generated the gif with the solutions. The scaling is a great tip, thanks. I'll try and faithfully reproduce your code for my stuff. Thanks a lot!!!!!! –  Amatya Sep 10 '12 at 20:21
can you please tell me how to extend this command eqs = (# == 0) & /@ (rhs/scale); to the case where rhs/scale is a 2 dimensional list. Thanks. If rhs/scale = {{y,y1},{z,z1}} then right now I am getting {{y,y1}==0,{z,z1}==0} but I would like to get {{y==0,y1==0},{z==0,z1==0}}. Thanks. –  Amatya Oct 2 '12 at 21:07