1
$\begingroup$

I´m new in coding and I need to solve a set of non linear equations. It´s a set of 3J equations in 2J + 2 unknowns: the J unknown values of βj, the J unknown values of μj, and the 2 unknown independent components of a unit vector u. I don´t know any reasonable starting values, so when I run my code I always got an error message that a single jacobian was encountered. The code is for J=2 but also I want to expand the system to J=10 (data already contained in variable mij.). Can anybody please help me solve these equations?

mij = {{1.3067138428220,2.4075852527272,3.338703707857},{1.1957072133028,1.1939672521530,1.6554978836236},{1.1634791565792,0.8419988047084135,1.167334669041159}, {1.1496251650317197,0.685390009666267,0.9499801680166198}, {1.1490966939110652,0.6854848022568772,0.95028283634603}, {1.1493538354209782,0.6903291660034211,0.9570670177527794}, {1.1379451243779293,0.568060100622182,0.7875227640478683},{1.132778443580569,0.5101431661797775,0.7071428683809919}, {1.1287596736205838,0.46716991388721535,0.6475558488667996}, {1.1218181817546011, 0.3911657652798617,0.5421184173117446}}; 
ni = {1.658747300215982695,5.231101511879049546,7.03671706263498898};
pi = {-0.041885130546511,-0.03337884298016919755,-0.016552882029393285823516};
FindRoot[
{ni[[1]] + muj1 ui1 - mij[[1, 1]] E^(-betaj1 pi[[1]]) == 0, 
ni[[2]] + muj1 ui1 - mij[[1, 2]] E^(-betaj1 pi[[2]]) == 0, 
ni[[3]] + muj1 ui1 - mij[[1, 3]] E^(-betaj1 pi[[3]]) == 0, 
ni[[1]] + muj2 ui2 - mij[[2, 1]] E^(-betaj2 pi[[1]]) == 0, 
ni[[2]] + muj2 ui2 - mij[[2, 2]] E^(-betaj2 pi[[2]]) == 0, 
ni[[3]] + muj2 ui2 - mij[[2, 3]] E^(-betaj2 pi[[3]]) == 0}, 
{{muj1, RandomReal[]}, {muj2, RandomReal[]}, {ui1, RandomReal[]}, {ui2,   RandomReal[]}, {betaj1, RandomReal[]}, {betaj2, RandomReal[]}}];

These are the original eq from the publication: original eq

Thanks for all your efforts, but some answers a kind of too sophisticated for me, I´m just a geologist ;) I´ve tried a lot, but still I cannot solve them. I believe, I don´t know how to arrange the eqs. In addition, the publication says that I can calculate the "2 unknown independent components" of the unit vector ui. When i=3, how do I get this "known component"? ui is supposed to be a unit vector along a mixing line, from which one side is known. I´ve tried a set up with j=3 to get 9 eqs and 9 variables, also got some more reasonable starting parameters, but I´m still stucked.

mij = {{1.3067138428220, 2.4075852527272, 
3.338703707857}, {1.1957072133028, 1.1939672521530, 
1.6554978836236}, {1.1634791565792, 0.8419988047084135, 
1.167334669041159}, {1.1496251650317197, 0.685390009666267, 
0.9499801680166198}, {1.1490966939110652, 0.6854848022568772, 
0.95028283634603}, {1.1493538354209782, 0.6903291660034211, 
0.9570670177527794}, {1.1379451243779293, 0.568060100622182, 
0.7875227640478683}, {1.132778443580569, 0.5101431661797775, 
0.7071428683809919}, {1.1287596736205838, 0.46716991388721535, 
0.6475558488667996}, {1.1218181817546011, 0.3911657652798617, 
0.5421184173117446}};
ni = {1.658747300215982695, 5.231101511879049546, 7.03671706263498898};
pi = {-0.041885130546511, -0.03337884298016919755, \
-0.016552882029393285823516};
FindRoot[
{
ni[[1]] + müj1 ui1 - mij2[[1, 1]] E^(-betaj1 pi[[1]]) == 0
, ni[[2]] + müj1 ui2 - mij2[[1, 2]] E^(-betaj1 pi[[2]]) == 0
, ni[[3]] + müj1 ui3 - mij2[[1, 3]] E^(-betaj1 pi[[3]]) == 0

, ni[[1]] + müj2 ui1 - mij2[[2, 1]] E^(-betaj2 pi[[1]]) == 0
, ni[[2]] + müj2 ui2 - mij2[[2, 2]] E^(-betaj2 pi[[2]]) == 0
, ni[[3]] + müj2 ui3 - mij2[[2, 3]] E^(-betaj2 pi[[3]]) == 0

, ni[[1]] + müj3 ui1 - mij2[[3, 1]] E^(-betaj3 pi[[1]]) == 0
, ni[[2]] + müj3 ui2 - mij2[[3, 2]] E^(-betaj3 pi[[2]]) == 0
, ni[[3]] + müj3 ui3 - mij2[[3, 3]] E^(-betaj3 pi[[3]]) == 0
}
, {{müj1, -0.150093}, {ui1, 0.48}, {ui2, 17.7990}, {betaj1, 
1.8125454}, {müj2, -0.2703640493}, {müj3, -0.27055}, {ui3, 
23.942731}, {betaj2, 1.801259006161}, {betaj3, 1.80586513960}}
, MaxIterations -> 10000]
$\endgroup$
  • 2
    $\begingroup$ You don't have enough unknowns as currently stated. Yes, you have six unknowns for six equations, but muj1 ui1 always appears in that combination, as does muj2 ui2. So you can replace those each with a single parameter. $\endgroup$ – KraZug Aug 25 '17 at 12:00
  • $\begingroup$ Thank, but sorry, I don´t understand it...when I replace the parameters I´ve got more eqs than unknowns plus I want to know muj and ui and not their products. I´ve posted the original eqs, maybe my coding of the eqs is the problem? $\endgroup$ – Stefan Aug 25 '17 at 12:40
  • $\begingroup$ Yes, your code is wrong. Each of your lines uses muj1 ui1, when you want those to change. $\endgroup$ – KraZug Aug 25 '17 at 12:47
  • 1
    $\begingroup$ Even if you fix your code, you still can't solve that equation as it stands - you need to decide with subset of equations you want to solve, such that #equations = #unknowns. $\endgroup$ – Lukas Lang Aug 25 '17 at 13:48
  • 1
    $\begingroup$ There is something amiss in that, as noted by @KraZug, it is not possible to separate out the factors in those mu-times-u products. So rechecking the equations is in order. Also, since there seem to be many more rows to mij, probably what is wanted is a best-fit rather than an exact solution to a subset (granted, the latter might be useful for initial guesses for the former). $\endgroup$ – Daniel Lichtblau Aug 25 '17 at 15:05
3
$\begingroup$

I will not re-copy the data nor the equations.

The second way you have written your equations is correct. You have 9 equations and 9 unknowns.

FindRoot is having trouble finding a solution.

Step 1 - FindMinimum

Recast the problem and use FindMinimum. We will take the left hand side of the equations, square them and add them, and then ask FindMinimum to determine the parameters that will minimize that expression.

This is more or less equivalent to setting each sub-expresion to zero.

solMin = FindMinimum[
  (ni[[1]] + müj1 ui1 - 
      mij2[[1, 1]] E^(-betaj1 pi[[1]]))^2 + (ni[[2]] + müj1 ui2 - 
      mij2[[1, 2]] E^(-betaj1 pi[[2]]))^2 + (ni[[3]] + müj1 ui3 - 
      mij2[[1, 3]] E^(-betaj1 pi[[3]]))^2 + (ni[[1]] + müj2 ui1 - 
      mij2[[2, 1]] E^(-betaj2 pi[[1]]))^2 + (ni[[2]] + müj2 ui2 - 
      mij2[[2, 2]] E^(-betaj2 pi[[2]]))^2 + (ni[[3]] + müj2 ui3 - 
      mij2[[2, 3]] E^(-betaj2 pi[[3]]))^2 + (ni[[1]] + müj3 ui1 - 
      mij2[[3, 1]] E^(-betaj3 pi[[1]]))^2 + (ni[[2]] + müj3 ui2 - 
      mij2[[3, 2]] E^(-betaj3 pi[[2]]))^2 + (ni[[3]] + müj3 ui3 - 
      mij2[[3, 3]] E^(-betaj3 pi[[3]]))^2,
  {
   {müj1, -0.150093}, {müj2, -0.2703640493}, {müj3, -0.27055},
   {ui1, 0.48}, {ui2, 17.7990}, {ui3, 23.942731},
   {betaj1, 1.8125454}, {betaj2, 1.801259006161}, {betaj3, 
    1.80586513960}
   }
  ]

(* {1.07279*10^-9, {müj1 -> -0.145231, müj2 -> -0.215279, 
  müj3 -> -0.23559, ui1 -> 1.71491, ui2 -> 18.4086, ui3 -> 24.7634, 
  betaj1 -> 1.81089, betaj2 -> 1.80424, betaj3 -> 1.80261}} *)

In order to validate the solution substitue the answer into the left hand side of the equations.

leftHandSide = {
   ni[[1]] + müj1 ui1 - mij2[[1, 1]] E^(-betaj1 pi[[1]]), 
   ni[[2]] + müj1 ui2 - mij2[[1, 2]] E^(-betaj1 pi[[2]]), 
   ni[[3]] + müj1 ui3 - mij2[[1, 3]] E^(-betaj1 pi[[3]]), 
   ni[[1]] + müj2 ui1 - mij2[[2, 1]] E^(-betaj2 pi[[1]]), 
   ni[[2]] + müj2 ui2 - mij2[[2, 2]] E^(-betaj2 pi[[2]]), 
   ni[[3]] + müj2 ui3 - mij2[[2, 3]] E^(-betaj2 pi[[3]]), 
   ni[[1]] + müj3 ui1 - mij2[[3, 1]] E^(-betaj3 pi[[1]]), 
   ni[[2]] + müj3 ui2 - mij2[[3, 2]] E^(-betaj3 pi[[2]]), 
   ni[[3]] + müj3 ui3 - mij2[[3, 3]] E^(-betaj3 pi[[3]])
   };

leftHandSide /. solMin[[2]]

(* {5.47011*10^-6, -7.00355*10^-6, 
 4.8275*10^-6, -8.14193*10^-6, 0.0000198635, -0.0000142023, 
 4.0679*10^-6, -0.0000138336, 0.000010002} *)

The answer is zero to within about 2*10-5.

Step 2 - Cross Correlation and FindRoot

There is a terrible cross-correlation between the parameters.

You can validate this by changing the starting points. You get a different answer but are still minimize the equations.

Here we will take the original FindRoot parameter and impose the starting value for müj1 in the first two equations.

solRoot = FindRoot[
  {
   ni[[1]] + -0.150093 ui1 - mij2[[1, 1]] E^(-betaj1 pi[[1]]) == 0, 
   ni[[2]] + -0.150093 ui2 - mij2[[1, 2]] E^(-betaj1 pi[[2]]) == 0, 
   ni[[3]] + müj1 ui3 - mij2[[1, 3]] E^(-betaj1 pi[[3]]) == 0, 
   ni[[1]] + müj2 ui1 - mij2[[2, 1]] E^(-betaj2 pi[[1]]) == 0, 
   ni[[2]] + müj2 ui2 - mij2[[2, 2]] E^(-betaj2 pi[[2]]) == 0, 
   ni[[3]] + müj2 ui3 - mij2[[2, 3]] E^(-betaj2 pi[[3]]) == 0, 
   ni[[1]] + müj3 ui1 - mij2[[3, 1]] E^(-betaj3 pi[[1]]) == 0, 
   ni[[2]] + müj3 ui2 - mij2[[3, 2]] E^(-betaj3 pi[[2]]) == 0, 
   ni[[3]] + müj3 ui3 - mij2[[3, 3]] E^(-betaj3 pi[[3]]) == 0
   },
  {
   {müj1, -0.150093}, {müj2, -0.2703640493}, {müj3, -0.27055},
   {ui1, 0.48}, {ui2, 17.7990}, {ui3, 23.942731},
   {betaj1, 1.8125454}, {betaj2, 1.801259006161}, {betaj3, 
    1.80586513960}
   },
  MaxIterations -> 10000
  ]

(* {müj1 -> -0.150102, müj2 -> -0.22251, müj3 -> -0.243507, 
 ui1 -> 1.65702, ui2 -> 17.8089, ui3 -> 23.9575, betaj1 -> 1.81694, 
 betaj2 -> 1.81301, betaj3 -> 1.81258} *)

Now test the result.

leftHandSide /. solRoot

(* {-0.0000153686, -0.000165174, 0., -5.55112*10^-17, \
-4.44089*10^-16, 0., 3.33067*10^-16, 0., 0.} *)

You could equally well have used -0.145231 for the value of müj1 and get an equally good answer.

solRoot = FindRoot[
  {
   ni[[1]] + -0.145231 ui1 - mij2[[1, 1]] E^(-betaj1 pi[[1]]) == 0, 
   ni[[2]] + -0.145231 ui2 - mij2[[1, 2]] E^(-betaj1 pi[[2]]) == 0, 
   ni[[3]] + müj1 ui3 - mij2[[1, 3]] E^(-betaj1 pi[[3]]) == 0, 
   ni[[1]] + müj2 ui1 - mij2[[2, 1]] E^(-betaj2 pi[[1]]) == 0, 
   ni[[2]] + müj2 ui2 - mij2[[2, 2]] E^(-betaj2 pi[[2]]) == 0, 
   ni[[3]] + müj2 ui3 - mij2[[2, 3]] E^(-betaj2 pi[[3]]) == 0, 
   ni[[1]] + müj3 ui1 - mij2[[3, 1]] E^(-betaj3 pi[[1]]) == 0, 
   ni[[2]] + müj3 ui2 - mij2[[3, 2]] E^(-betaj3 pi[[2]]) == 0, 
   ni[[3]] + müj3 ui3 - mij2[[3, 3]] E^(-betaj3 pi[[3]]) == 0
   },
  {
   {müj1, -0.145231}, {müj2, -0.2703640493}, {müj3, -0.27055},
   {ui1, 0.48}, {ui2, 17.7990}, {ui3, 23.942731},
   {betaj1, 1.8125454}, {betaj2, 1.801259006161}, {betaj3, 
    1.80586513960}
   },
  MaxIterations -> 10000
  ]

(* {müj1 -> -0.14524, müj2 -> -0.215302, müj3 -> -0.235619, 
 ui1 -> 1.71249, ui2 -> 18.4051, ui3 -> 24.7595, betaj1 -> 1.81694, 
 betaj2 -> 1.81301, betaj3 -> 1.81258} *)

Virtually identical results with a different müj1.

Summary

The problem is ill-conditioned. You might as well impose the value of one of the product terms and then use FindRoot to determine the other values.

There is not sufficient information to extract nine independent parameters from the data and equations.

$\endgroup$
  • $\begingroup$ Hi Guys, hi Jack! Thanks for all your efforts. I´ve processed the data from Jacks solution through the upcoming steps, but the results are desperate incorrect. I will try to find another approach for my work than solving these horrible eqs. Thanks @all! $\endgroup$ – Stefan Aug 28 '17 at 11:11

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.