I recently tried to find a root with FindRoot
in Mathematica. I know that the problem has a solution, but when using FindRoot, this error appears :
FindRoot::lstol: The line search decreased the step size to within tolerance specified by AccuracyGoal and PrecisionGoal but was unable to find a sufficient decrease in the merit function. You may need more than 10000.` digits of working precision to meet these tolerances. >> t->{0.1139}
However, I get an answer, but this is not a really accurate one! Even though I changed AccuracyGoal
, PrecessionGoal
and WorkingPrecession
, the answer didn't change. Please help me getting rid of the error.
My code is:
a1 = 0.375;
a2 = 0.175;
b1 = 0.125;
b2 = 0.325;
m1 = 0.05;
m2 = 0.5;
mh = 0.5;
l1 = a1 + b1;
l2 = a2 + b2;
λ = 0.0504;
g = 9.81;
kinetic1 =
1/2 m1 (a1^2 Cos[q1[t]]^2 q1'[t]^2 + a1^2 Sin[q1[t]]^2 q1'[t]^2) +
1/2 m2 ((a2 Cos[q1[t]] q1'[t] +
l1 Cos[q1[t]] q1'[t])^2 + (-a2 Sin[q1[t]] q1'[t] -
l1 Sin[q1[t]] q1'[t])^2) +
1/2 mh ((l1 Cos[q1[t]] q1'[t] +
l2 Cos[q1[t]] q1'[t])^2 + (-l1 Sin[q1[t]] q1'[t] -
l2 Sin[q1[t]] q1'[t])^2) +
1/2 m2 ((l1 Cos[q1[t]] q1'[t] + l2 Cos[q1[t]] q1'[t] +
b2 Cos[q2[t]] q2'[t])^2 + (-l1 Sin[q1[t]] q1'[t] -
l2 Sin[q1[t]] q1'[t] - b2 Sin[q2[t]] q2'[t])^2) +
1/2 m1 ((l1 Cos[q1[t]] q1'[t] + l2 Cos[q1[t]] q1'[t] +
l2 Cos[q2[t]] q2'[t] +
b1 Cos[q3[t]] q3'[t])^2 + (-l1 Sin[q1[t]] q1'[t] -
l2 Sin[q1[t]] q1'[t] - l2 Sin[q2[t]] q2'[t] -
b1 Sin[q3[t]] q3'[t])^2);
potential1 =
m1 g (a1 Cos[q1[t]]) + m2 g ((l1 + a2) Cos[q1[t]]) +
mh g ((l1 + l2) Cos[q1[t]]) +
m2 g ((l1 + l2) Cos[q1[t]] - b2 Cos[q2[t]]) +
m1 g ((l1 + l2) Cos[q1[t]] - l2 Cos[q2[t]] - b1 Cos[q3[t]]);
lagrangeq11 =
D[D[kinetic1, q1'[t]], t] - D[kinetic1, q1[t]] +
D[potential1, q1[t]];
lagrangeq21 =
D[D[kinetic1, q2'[t]], t] - D[kinetic1, q2[t]] +
D[potential1, q2[t]];
lagrangeq31 =
D[D[kinetic1, q3'[t]], t] - D[kinetic1, q3[t]] +
D[potential1, q3[t]];
kinetic2 =
1/2 m1 (a1^2 Cos[q1[t]]^2 q1'[t]^2 + a1^2 Sin[q1[t]]^2 q1'[t]^2) +
1/2 m2 ((a2 Cos[q1[t]] q1'[t] +
l1 Cos[q1[t]] q1'[t])^2 + (-a2 Sin[q1[t]] q1'[t] -
l1 Sin[q1[t]] q1'[t])^2) +
1/2 mh ((l1 Cos[q1[t]] q1'[t] +
l2 Cos[q1[t]] q1'[t])^2 + (-l1 Sin[q1[t]] q1'[t] -
l2 Sin[q1[t]] q1'[t])^2) +
1/2 m2 ((l1 Cos[q1[t]] q1'[t] + l2 Cos[q1[t]] q1'[t] +
b2 Cos[q2[t]] q2'[t])^2 + (-l1 Sin[q1[t]] q1'[t] -
l2 Sin[q1[t]] q1'[t] - b2 Sin[q2[t]] q2'[t])^2) +
1/2 m1 ((l1 Cos[q1[t]] q1'[t] + l2 Cos[q1[t]] q1'[t] +
l2 Cos[q2[t]] q2'[t] +
b1 Cos[q2[t]] q2'[t])^2 + (-l1 Sin[q1[t]] q1'[t] -
l2 Sin[q1[t]] q1'[t] - l2 Sin[q2[t]] q2'[t] -
b1 Sin[q2[t]] q2'[t])^2);
potential1 =
m1 g (a1 Cos[q1[t]]) + m2 g ((l1 + a2) Cos[q1[t]]) +
mh g ((l1 + l2) Cos[q1[t]]) +
m2 g ((l1 + l2) Cos[q1[t]] - b2 Cos[q2[t]]) +
m1 g ((l1 + l2) Cos[q1[t]] - l2 Cos[q2[t]] - b1 Cos[q2[t]]);
lagrangeq12 =
D[D[kinetic2, q1'[t]], t] - D[kinetic2, q1[t]] +
D[potential2, q1[t]];
lagrangeq22 =
D[D[kinetic2, q2'[t]], t] - D[kinetic2, q2[t]] +
D[potential2, q2[t]];
q10 = -0.1877 - λ;
q30 = q20 = 0.2884 - λ;
q1d0 = 1.1014;
q3d0 = q2d0 = 0.0399;
beforeknee1 = {{0, q10}}; beforeknee2 = {{0, q20}}; beforeknee3 = {{0,
q30}};
time = timep = 0;
Clear[v1];
Clear[v2];
lagrangesolveg1 =
NDSolve[{lagrangeq11 == 0, lagrangeq21 == 0, lagrangeq31 == 0,
q1[0] == q10, q2[0] == q20, q3[0] == q20, q1'[0] == q1d0,
q2'[0] == q2d0, q3'[0] == q2d0}, {q1, q3, q2}, {t, 0, 1},
MaxSteps -> 1000000];
time1 = FindRoot[(q3[t] == q2[t]) /. lagrangesolveg1, {t, 0.5, 0,
1}]; time1 = t /. time1; time1 = SetAccuracy[time1, 5];
For[τ = 1, τ < (IntegerPart[10000 time1] + 1), τ++,
q1t = q1[τ/10000] /. lagrangesolveg1;
q2t = q2[τ/10000] /. lagrangesolveg1;
q3t = q3[τ/10000] /. lagrangesolveg1;
time = N[timep + τ/10000];
beforeknee1 = Append[beforeknee1, {time, q1t[[1]]}];
beforeknee2 = Append[beforeknee2, {time, q2t[[1]]}];
beforeknee3 = Append[beforeknee3, {time, q3t[[1]]}];];
timep = time;
hh1 = q1'[time1] /. lagrangesolveg1;
hh1 = hh1[[1]];
hh2 = q2'[time1] /. lagrangesolveg1;
hh2 = hh2[[1]];
hh3 = q3'[time1] /. lagrangesolveg1;
hh3 = hh3[[1]];
hh4 = q1[time1] /. lagrangesolveg1;
hh4 = hh4[[1]];
hh5 = q3[time1] /. lagrangesolveg1;
hh5 = hh5[[1]];
kneeright1 =
Cross[{a1 Sin[hh4], a1 Cos[hh4], 0},
m1 {a1 Cos[hh4] v1, -a1 Sin[hh4] v1, 0}] +
Cross[{(a2 + l1) Sin[hh4], (a2 + l1) Cos[hh4], 0},
m2 {a2 Cos[hh4] v1 + l1 Cos[hh4] v1, -a2 Sin[hh4] v1 -
l1 Sin[hh4] v1, 0}] +
Cross[{(l2 + l1) Sin[hh4], (l2 + l1) Cos[hh4], 0},
mh {l1 Cos[hh4] v1 + l2 Cos[hh4] v1, -l1 Sin[hh4] v1 -
l2 Sin[hh4] v1, 0}] +
Cross[{(l2 + l1) Sin[hh4] - b2 Sin[hh5], (l2 + l1) Cos[hh4] -
b2 Cos[hh5], 0},
m2 {l1 Cos[hh4] v1 + l2 Cos[hh4] v1 +
b2 Cos[hh5] v2, -l1 Sin[hh4] v1 - l2 Sin[hh4] v1 -
b2 Sin[hh5] v2, 0}] +
Cross[{(l2 + l1) Sin[hh4] - l2 Sin[hh5] -
b1 Sin[hh5], (l2 + l1) Cos[hh4] - l2 Cos[hh5] - b1 Cos[hh5], 0},
m1 {l1 Cos[hh4] v1 + l2 Cos[hh4] v1 + l2 Cos[hh5] v2 +
b1 Cos[hh5] v2, -l1 Sin[hh4] v1 - l2 Sin[hh4] v1 -
l2 Sin[hh5] v2 - b1 Sin[hh5] v2, 0}];
kneeleftt1 =
Cross[{a1 Sin[hh4], a1 Cos[hh4], 0},
m1 {a1 Cos[hh4] hh1, -a1 Sin[hh4] hh1, 0}] +
Cross[{(a2 + l1) Sin[hh4], (a2 + l1) Cos[hh4], 0},
m2 {a2 Cos[hh4] hh1 + l1 Cos[hh4] hh1, -a2 Sin[hh4] hh1 -
l1 Sin[hh4] hh4, 0}] +
Cross[{(l2 + l1) Sin[hh4], (l2 + l1) Cos[hh4], 0},
mh {l1 Cos[hh4] hh1 + l2 Cos[hh4] hh1, -l1 Sin[hh4] hh1 -
l2 Sin[hh4] hh1, 0}] +
Cross[{(l2 + l1) Sin[hh4] - b2 Sin[hh5], (l2 + l1) Cos[hh4] -
b2 Cos[hh5], 0},
m2 {l1 Cos[hh4] hh1 + l2 Cos[hh4] hh1 +
b2 Cos[hh5] hh2, -l1 Sin[hh4] hh1 - l2 Sin[hh4] hh1 -
b2 Sin[hh5] hh2, 0}] +
Cross[{(l2 + l1) Sin[hh4] - l2 Sin[hh5] -
b1 Sin[hh5], (l2 + l1) Cos[hh4] - l2 Cos[hh5] - b1 Cos[hh5], 0},
m1 {l1 Cos[hh4] hh1 + l2 Cos[hh4] hh1 + l2 Cos[hh5] hh2 +
b1 Cos[hh5] hh3, -l1 Sin[hh4] hh1 - l2 Sin[hh4] hh1 -
l2 Sin[hh5] hh2 - b1 Sin[hh2] hh3, 0}];
kneeright2 =
Cross[{-b2 Sin[hh5], -b2 Cos[hh5], 0},
m2 {l1 Cos[hh4] v1 + l2 Cos[hh4] v1 +
b2 Cos[hh5] v2, -l1 Sin[hh4] v1 - l2 Sin[hh4] v1 -
b2 Sin[hh5] v2, 0}] +
Cross[{-l2 Sin[hh5] - b1 Sin[hh5], -l2 Cos[hh5] - b1 Cos[hh5], 0},
m1 {l1 Cos[hh4] v1 + l2 Cos[hh4] v1 + l2 Cos[hh5] v2 +
b1 Cos[hh5] v2, -l1 Sin[hh4] v1 - l2 Sin[hh4] v1 -
l2 Sin[hh5] v2 - b1 Sin[hh5] v2, 0}];
kneeleftt2 =
Cross[{-b2 Sin[hh5], -b2 Cos[hh5], 0},
m2 {l1 Cos[hh4] hh1 + l2 Cos[hh4] hh1 +
b2 Cos[hh5] hh2, -l1 Sin[hh4] hh1 - l2 Sin[hh4] hh1 -
b2 Sin[hh5] hh2, 0}] +
Cross[{-l2 Sin[hh5] - b1 Sin[hh5], -l2 Cos[hh5] - b1 Cos[hh5], 0},
m1 {l1 Cos[hh4] hh1 + l2 Cos[hh4] hh1 + l2 Cos[hh5] hh2 +
b1 Cos[hh5] hh3, -l1 Sin[hh4] hh1 - l2 Sin[hh4] hh1 -
l2 Sin[hh5] hh2 - b1 Sin[hh2] hh3, 0}];
kneestrike =
Solve[{kneeleftt1 == kneeright1, kneeleftt2 == kneeright2}, {v1,
v2}];
{{v1, v2}} = {v1, v2} /. kneestrike;
lagrangesolveg2 =
NDSolve[{lagrangeq12 == 0, lagrangeq22 == 0, q1[0] == hh4,
q2[0] == hh5, q1'[0] == v1, q2'[0] == v2}, {q1, q2}, {t, 0, 1},
MaxSteps -> 1000000];
time2 = FindRoot[(q1[t] == -q2[t]) /. lagrangesolveg2, {t, 0.5, 0, 1}]