How to transform a FindRoot in a NDSolve?

Question

some days ago I did a question about a problem with a findroot. The solution was to convert the FindRoot in a "fake" NDSolve. How can I do the same if the FindRoot has more equation ? I make here an example, this is the FindRoot version:

tabexp=ParallelTable[{y1,y2/.FindRoot[Pinj-Fun[y1,y2],{y2,0.2,2},AccuracyGoal->25,PrecisionGoal->25,MaxIterations->2000]},{y1,0.2,10,0.9/200}];

This the NDSolve versione:

 sol=NDSolve[{Pinj==Fun[y1,y2[y1]],x'[y1]==1,x[0.2]==0.2,y2[0.2]==(y2y2/.FindRoot[Pinj-Fun[0.2,y2y2],{y2y2,2}])},y2,{y1,0.2,10}]

This is more or less clear to me, the problem is that now i would like to solve in this way a FindRoot with two equations:

 tabexp=ParallelTable[{FindRoot[{Pinjx-Fun[y1,y2]- Loss[y1,1],Fun[y1,y2]-Loss[y2,1]},{y1,5},{y2,5},AccuracyGoal->25,PrecisionGoal->25,MaxIterations->2000]},{Pinjx,0.2,10,9.8/50}];

The rest of the code is here

 n[ω_?NumericQ,y_?NumericQ]:=(Exp[ ω/y]-1)-1;
Fun[y1_?NumericQ,y2_?NumericQ]:= NIntegrate[ω (n[ω,y1]-n[ω,y2]),{ω,0,10},MinRecursion->4,Method->{"GlobalAdaptive","MaxErrorIncreases"->100000,"SymbolicProcessing"->0,"SingularityHandler"->None},PrecisionGoal->5,WorkingPrecision->16];
Loss[y_,ye_]:=y5-ye5;
Pinj=0.01;

score 2 · Accepted Answer · answered Nov 21 '18 at 02:43

Solving multiple equations as a DAE is similar to solving just one. Put them all in:

n[ω_?NumericQ, y_?NumericQ] := (Exp[ω/y] - 1) - 1;  (* sure this is right? *)
Fun[y1_?NumericQ, y2_?NumericQ] := 
  NIntegrate[ω (n[ω, y1] - n[ω, y2]), {ω, 
    0, 10}, Method -> {"GaussKronrodRule", "Points" -> 11}, 
   PrecisionGoal -> 10];
Loss[y_, ye_] := y - ye;
Pinj = 1/100;

Block[{Pinjx = 2/10},
 yy0 = {y1, y2} /.   (* initial starting values for y1, y2 *)
   FindRoot[{Pinjx - Fun[y1, y2] == Loss[y1, 1], 
     Fun[y1, y2] == Loss[y2, 1]}, {{y1, 2}, {y2, 3}}]
 ]

{ndsol} = NDSolve[
  {Pinjx - Fun[y1[Pinjx], y2[Pinjx]] == Loss[y1[Pinjx], 1], 
   Fun[y1[Pinjx], y2[Pinjx]] == Loss[y2[Pinjx], 1],
   {y1[2/10], y2[2/10]} == yy0,
   x'[Pinjx] == 1, x[2/10] == 2/10},
  {y1, y2}, {Pinjx, 2/10, 10}]

(I tweaked the NIntegrate code to make it converge more reliably & faster.)

How to transform a FindRoot in a NDSolve?

1 Answers1

Linked