Solving a system of 5 dependent non-linear equations

I am trying to solve the system of 5 non-linear equations:

Io * exp(Voc/ (m*Vt)) + Voc/Rsh - Is = 0

(Isc - (Voc - Rs*Isc)/ Rsh) * exp(-Voc/ (m*Vt)) - Io = 0

Isc - (Vmp + Rs*Imp + Rs*Isc)/ Rsh - (Isc - (Voc - Rs*Isc)/ Rsh) * exp((Vmp + Rs*Imp - Voc)/ (m*Vt)) - Imp = 0

Imp + ((-(Rsh*Isc - Voc + Rs*Isc) * exp((Vmp + Rs*Imp - Voc)/ (m*Vt))/ (Rsh * m*Vt) - 1/Rsh)/ (1 + Rs * (Rsh*Isc - Voc + Rs*Isc) * exp((Vmp + Rs*Imp - Voc)/ (m*Vt))/ (Rsh * m*Vt) + Rs/Rsh)) * Vmp = 0

(-(Rsh*Isc - Voc + Rs*Isc) * exp((Rs*Isc - Voc)/ (m*Vt))/ (Rsh * m*Vt) - 1/Rsh)/ (1 + Rs * (Rsh*Isc - Voc + Rs*Isc) * exp((Rs*Isc - Voc)/ (m*Vt))/ (Rsh * m*Vt) + Rs/Rsh) + 1/Rsh = 0

in which m, Rsh, Rs, Io and Is are the unknown variables. The other variables have the values:

Vmp = 31.1;
Imp = 8.67;
Voc = 38.2;
Isc = 9.19;
K = 1.38 * 10^-23;
T = 298.15;
q = 1.6 * 10*-19; 
Vt = K*T/q;

I know the solution must be around m=62.3, Rsh=786, Rs=0.2748, Io=407.33, Is=9.1932

I have tried using the function fsolve as shown below

syms  m Rsh Rs Io Is
V = [m, Rsh, Rs, Io, Is];
 F = @(V) [V(4) * exp(Voc/ (V(1)*Vt)) + Voc/V(2) - V(5);
          (Isc - (Voc - V(3)*Isc)/ V(2)) * exp(-Voc/ (V(1)*Vt)) - V(4);
          Isc - (Vmp + V(3)*Imp + V(3)*Isc)/ V(2) - (Isc - (Voc - V(3)*Isc)/ V(2)) * exp((Vmp + V(3)*Imp - Voc)/ (V(1)*Vt)) - Imp;
          Imp + ((-(V(2)*Isc - Voc + V(3)*Isc) * exp((Vmp + V(3)*Imp - Voc)/ (V(1)*Vt))/ (V(2) * V(1)*Vt) - 1/V(2))/ (1 + V(3) * (V(2)*Isc - Voc + V(3)*Isc) * exp((Vmp + V(3)*Imp - Voc)/ (V(1)*Vt))/ (V(2) * V(1)*Vt) + V(3)/V(2))) * Vmp;
          (-(V(2)*Isc - Voc + V(3)*Isc) * exp((V(3)*Isc - Voc)/ (V(1)*Vt))/ (V(2) * V(1)*Vt) - 1/V(2))/ (1 + V(3) * (V(2)*Isc - Voc + V(3)*Isc) * exp((V(3)*Isc - Voc)/ (V(1)*Vt))/ (V(2) * V(1)*Vt) + V(3)/V(2)) + 1/V(2)];
InitialGuess = [50; 100; 1; 100; 10];
XY = fsolve(F, InitialGuess);

but I get the Error using trustnleqn (line 28). Objective function is returning undefined values at initial point. FSOLVE cannot continue. error message. Can someone tell how they would solve this problem?

[John D'Errico]

As I said in my answer, fsolve will NEVER be able to solve this. Even if you learn to use fsolve properly it will fail.

So, what are you doing wrong with fsolve?

DON'T USE SYMS!!!!!!!!!!

People think that in order to make a function that has an unknown in it, they need to use symbolic variables. NOT TRUE!

For example, the function F below has an unknown value in it, thus x. I never made x symbolic. There is no need to do so. You can use a root finder on F, for example, fzero.

F = @(x) x.^2 - 2;
fzero(F,1)
ans =
     1.4142

Regardless, you might try reading the examples for fsolve. Read the docs.

doc fsolve

Anyway, fsolve is irrelevant. You CANNOT solve this using fsolve, so spending the time teaching you what you did wrong is a complete waste of time here. Sorry, but it is. If you want to learn how to use fsolve on your own time, again, try the examples provided for fsolve.

As I said in my answer, the problem is so poorly conditioned that you CANNOT solve it using double precision arithmetic. In fact, even trying to solve it using tools like vpasolve will fail, if you are correct about the value of m.

What you seem not to recognize is that a number like

exp(-1e22)

is unimaginably small. How many zeros does it have before you get to the first non-zero digit? Not just 22 zero. You cannot write all of those zeros in your life time. Your computer cannot even count them in your life time.

[Frederico Lourenço]

I'm sorry, I understand your frustration, I really don't know much about how most of the functions work, that's why I came here. Maybe someone would be able to explain to me what I was doing wrong.

Anyway, that unimaginably small number that you refer was due to a typo in the value of q. The real value for Vt is 0.0257, and not -1.3534e-23, which I thank you for pointing to me. I already corrected that and since then I no longer get the undefined values that crash fsolve.

My problem now, like I mentioned in your answer, is that the function just does not converge and as you said I don't understand what I'm doing wrong since I'm following the explanation in the link of the MathWorks Support Team.

You say that I should not use fsolve. Well my original question was how would someone solve the problem I presented, solve a system of 5 dependent non-linear equations.

Answer 1

[John BG]

Hi Frederico

this is John BG <mailto:[e-mail] [e-mail]>

1.

Actually, command fsolve works pretty well here, let me explain:

Since V(5) unknown only shows up in the 1st equation I decided to put off the 1st equation, let's solve the remaining system with fsolve

clear all;clc;close all

syms m Rsh Rs Io Is

Vmp = 31.1;
Imp = 8.67;
Voc = 38.2;
Isc = 9.19;
K = 1.38 * 10^-23;
T = 298.15;
q = -1.6 * 10^-19; 
Vt = K*T/q;

V = [m, Rsh, Rs, Io];
 F = @(V) [  (Isc - (Voc - V(3)*Isc)/ V(2)) * exp(-Voc/ (V(1)*Vt)) - V(4);
          Isc - (Vmp + V(3)*Imp + V(3)*Isc)/ V(2) - (Isc - (Voc - V(3)*Isc)/ V(2)) * exp((Vmp + V(3)*Imp - Voc)/ (V(1)*Vt)) - Imp;
          Imp + ((-(V(2)*Isc - Voc + V(3)*Isc) * exp((Vmp + V(3)*Imp - Voc)/ (V(1)*Vt))/ (V(2) * V(1)*Vt) - 1/V(2))/ (1 + V(3) * (V(2)*Isc - Voc + V(3)*Isc) * exp((Vmp + V(3)*Imp - Voc)/ (V(1)*Vt))/ (V(2) * V(1)*Vt) + V(3)/V(2))) * Vmp;
          (-(V(2)*Isc - Voc + V(3)*Isc) * exp((V(3)*Isc - Voc)/ (V(1)*Vt))/ (V(2) * V(1)*Vt) - 1/V(2))/ (1 + V(3) * (V(2)*Isc - Voc + V(3)*Isc) * exp((V(3)*Isc - Voc)/ (V(1)*Vt))/ (V(2) * V(1)*Vt) + V(3)/V(2)) + 1/V(2)];
InitialGuess = [50; 0.1; .01; 1];

options = optimoptions('fsolve','Display','none','PlotFcn',@optimplotfirstorderopt,'MaxFunctionEvaluations',2000);
XY = fsolve(F, InitialGuess,options);

XY =
  51.755945628868020
   3.071518542695716
   1.085173513887406
  -1.062922141399835

2.

These are semiconductor equations, I like reading what is what, so I am putting back the initial variable names

m=XY(1)
Rsh=XY(2)
Rs=XY(3)
Io=XY(4)

3.

Now the initial equation is trivial

Is=  Io * exp(Voc/ (m*Vt) )+ Voc/Rsh

Is =
  12.436844990189343

4.

fsolve even produces the optimality graph ending up in 3.5e9 quite high.

5.

I was wondering, should your 1st equation have an inverted polarity, then Is would be

Is=  Io * exp(-Voc/ (m*Vt) )+ Voc/Rsh
Is =
    -3.101249033198816e+12

6.

Please note that function fsolve has a default amount if iterations equal to 400. If not changed, for this equations system, fsolve completes all the loops, without a solution, but warns that it may have a solution if increasing such default value.

To so, just change in advance the options field

'MaxFunctionEvaluations'

to for instance

2000

.

Frederico

Cosi fan tutte :)

If you find this answer useful would you please be so kind to consider marking my answer as Accepted Answer?

To any other reader, if you find this answer useful please consider clicking on the thumbs-up vote link

thanks in advance

John BG

<mailto:[e-mail] [e-mail]>

[Frederico Lourenço]

Thanks a lot for your answer, John! It was really helpful.

Since, as you noted in point 5, the system is not really convergent I am assuming that it is because there is something wrong with the equations, maybe a typo somewhere that is preventing the system from getting to the correct result. I will check once again if the equations are correct and will comment here if you solved my problem or not!

EDIT: The equations are correct. However, since the last 3 equations can be calculated separately to get m, Rs and Rsh, we can heavily simplify the problem to only 3 equations with 3 unknowns. This helps the system a lot, and it becomes convergent! The last 2 equations can be used afterwards, to calculate Io and then Is directly, and the results end up very close to what I was predicting.

I thank you guys a lot for your time and patience, John D'Errico and John BG!

[John BG]

At first sight, I don't see syntax error, but it may be misleading because the equations may not be the right ones.

You are modelling a semiconductor, what is it? a diode? a transistor? can you show the small signal model? with all capacitors, controlled sources and whatever you have used to get to the equations shown above?

Having a look at the transfer function mat tell if it's going to burn when you switch it on.

If you have S parameters, then check for K, the stability factor on the Smith chart for different frequencies. You will see whether the work point is in safe area or again there's oscillation, that may be the case. Are you designing an oscillator?

[Frederico Lourenço]

Sorry, I said stable and I meant convergent. It is an equivalent model of an electrical circuit with a diode and two resistors, no oscillatory behavior

Answer 2

[John D'Errico]

Maybe you THINK you know where the solution lies. But if you are correct, you will essentially NEVER be able to solve this problem using tools like fsolve.

I'll use symbolic tools here.

syms m Rsh Rs Io Is 
E(1) = Io * exp(Voc/ (m*Vt)) + Voc/Rsh - Is == 0;
E(2) = (Isc - (Voc - Rs*Isc)/ Rsh) * exp(-Voc/ (m*Vt)) - Io == 0;

Etc.

We can stop at equation 1 in fact.

pretty(vpa(E(1),5))
      /            24 \
      |   2.8224 10   |            38.2
Io exp| - ----------- | - 1.0 Is + ---- == 0.0
      \        m      /             Rsh

You claim that m is around 62.3.

So that would have you compute exp(-2.8224e24/62.3). The number is so small that it is effectively zero, even when using symbolic or my own HPF tools.

In fact, every place where you have m, you divide it into a number on the order of 1e23 or so, then you exponentiate. So m is irrelevant. It results in an exponential that is zero.

pretty(vpa(E(2),5))
   /          24 \
   | 2.8224 10   | / 9.19 Rs - 38.2        \
exp| ----------- | | -------------- + 9.19 | - 1.0 Io == 0.0
   \      m      / \       Rsh             /

>> pretty(vpa(E(3),5))
                                  /                            22 \
       (17.86 Rs + 31.1) 1.0      |   (8.67 Rs - 7.1) 7.3886 10   | / 9.19 Rs - 38.2        \
0.52 - --------------------- - exp| - --------------------------- | | -------------- + 9.19 | 1.0 == 0.0
                Rsh               \                m              / \       Rsh             /

So, even if you are calling fsolve wrongly, you are wasting your time. In fact, you are wasting your time trying to solve this even with symbolic tools, IF m is on the order of 60. Unless m was on the order of 1e20 or more, you will have no chance of solving this. And even then, fsolve would be a complete waste of time.

[Frederico Lourenço]

First of all, thanks for your fast answer! What you say is indeed true, I had a typo in the value of q. It was suppose to be 1.6 * 10^(-19), and not 1.6 * 10*-19.

The error went away when I changed that but unfortunately I believe I still have a mistake because when I do the line

ShouldBeZero = F(XY)

like in Is it possible to solve simultaneous nonlinear equations with MATLAB? (link), the results are

ShouldBeZero =
   -2.2483
    0.0165
   -3.0388
    0.8064
   -0.0109

The program only stops because of this

fsolve stopped because it exceeded the function evaluation limit,
options.MaxFunctionEvaluations = 500 (the default value).

and not because it has reached results within the tolerance of the function. Do you know why this happens?

Answer 3

[Alex Sha]

Hi, if vt = 0.0257, the solution should be:

m: 512.891368076595
rsh: -1.08505661410065
rs: 3.40678689852078
io: 0.856812866408915
is: -19.6641364483436

Fevl:
3.80140363631654E-13
1.59949831157746E-12
4.13891143580258E-13
-2.0961010704923E-13
3.41060513164848E-13

while, if Vt = K*T/q=-1.35344407894737E-23, the solution will be:

m: -9.73908591023589E23
rsh: -1.08505661410068
rs: 3.40678689852086
io: 0.856812866410473
is: -19.6641364483141

Fevl:
1.4210854715202E-14
1.11022302462516E-16
-8.34887714518118E-14
-1.24344978758018E-14
-1.55431223447522E-15

[Frederico Lourenço]

It's been two years since I posted this and got the answer, but thank you for your comment.

Answer 4

[Alex Sha]

My pleasure, it is actually very challenging problem, hard to be solved by common software or algorithms.