Take the 2-minute tour ×
Mathematica Stack Exchange is a question and answer site for users of Mathematica. It's 100% free, no registration required.

We are considering the following equations

m = M/2 (1 + Tanh[v[x]/vo]);
q = Q/2 (1 + Tanh[v[x]/vo]);
M = 2;
Q = 1;
vo = 0.01;
d = 4;
small = 1/10000;
Eqn1 = 1/2 d z[x]^d  D[v[x], x]^2  m - z[x]^d  D[v[x], x]^2  m - 
   d z[x]^(2  d - 2)  D[v[x], x]^2  q^2 + 
   2 z[x]^(2  d - 2)  D[v[x], x]^2  q^2 + z[x]  D[v[x], {x, 2}] + 
   2  D[v[x], x]  D[z[x], x] + D[v[x], x]^2 - 1;
Eqn2 = -(1/2)  z[x]^d  D[v[x], x]^2  D[m, v[x]] - 
   z[x]^d  D[v[x], {x, 2}]  m - 
   d z[x]^( d - 1)   D[v[x], x]  D[z[x], x]   m + 
   z[x]^(2  d - 2)   D[v[x], x]^2  q  D[ q, z[x]]   + 
   z[x]^(2  d - 2)    D[v[x], {x, 2}]   q^2 + 
   2  d  z[x]^(2  d - 3)  D[v[x], x]  D[z[x], x]   q^2 - 
   2   z[x]^(2  d - 3)    D[v[x], x]  D[z[x], x]  q^2 +  
   D[v[x], {x, 2}] + D[z[x], {x, 2}];

When we solve the two equations with a pair of zs, and vs, for example

zs = 0.85; vs = -0.04529;  

we can get the solutions though there is a singularity :

NDSolve::ndsz: At x == -1.00007, step size is effectively zero; singularity or stiff system suspected. >>

But when we run the program

l = 3; zsi = 0.81; zse = 1.5; vsi = -1.5; vse = 0.5; z0 = 0.01;
Table[sol = NDSolve[{Eqn1 == 0, Eqn2 == 0, z[0] == zs,
    v[0] == vs,
    z'[0] == 0,
    v'[0] == 0}, {z, v}, {x, -1, 1}, MaxStepFraction -> 1/11, 
   Method -> {"BDF", "MaxDifferenceOrder" -> 4}]; 
 evz[x_] := First[z[x] /. sol];
 k = FindRoot[evz[lh] - z0, {lh, zs}, MaxIterations -> 1000]; 
 long = 2*lh /. k; CCC = evz[long];
 If[l*0.9995 < long < 
   l*1.0005, {v[long/2] /. 
    sol, (2 NIntegrate[zs/(z[x] /. sol)^2, {x, 0, l/2}] - 
      2  Log[2/0.01])/l}, {0, 0}], {zs, zsi, zse, 0.05}, {vs, vsi, vse, 0.05}]

the singularity will stop the program to run. In this program, we Use the Table command to vary zs and vs, and want to choose a pair of zs and vs that satisfy the condition in the If command. So we must find a way to cancel out the singularity and run our program. This problem has perplexed me several months.


IN the First cell, i want to solve Eq1 and Eq2, namely

   sol = NDSolve[{Eqn1 == 0, Eqn2 == 0, z[0] == zs,
   v[0] == vs,
   z'[0] == 0,
   v'[0] == 0}, {z, v}, {x, -1.5, 1.5}, MaxStepFraction -> 1/11, 
  Method -> {"BDF", "MaxDifferenceOrder" -> 4}]

You can run our program, especially the second cell, you will find some numerical problems. I want to solve these problems to run the second cell

share|improve this question

closed as off-topic by Michael E2, rasher, Oleksandr R., m_goldberg, Yves Klett Apr 20 at 7:21

  • The question does not concern the technical computing software Mathematica by Wolfram Research. Please see the help center to find out about the topics that can be asked here.
If this question can be reworded to fit the rules in the help center, please edit the question.

1  
"If i have not narrated this problem clearly..." - then edit your question instead of asking to be e-mailed. –  J. M. Oct 5 '12 at 13:25
    
How did you solve the First cell in your question? Through FindRoot? And what are you trying to do in the second cell in your question? You seem to be using a stiff solver which should have eased your troubles... briefly.. –  drN Oct 6 '12 at 3:50
    
I think you can still simplify your question a little, for example, CCC = evz[long]; is useless in the sample, right? –  xzczd Oct 10 '12 at 11:08
4  
This question appears to be off-topic because it is about a singularity in the mathematical model and not about Mathematica; further the OP has been absent for over a year. –  Michael E2 Apr 19 at 19:19

1 Answer 1

This question appears to be off-topic because it is about a singularity in the mathematical model and not about Mathematica; further the OP has been absent for over a year.

For example:

l = 3; zsi = 0.81; zse = 1.5; vsi = -1.5; vse = 0.5; z0 = 0.01;
Block[{zs = zsi, vs = vsi},
  sol = NDSolve[{Eqn1 == 0, Eqn2 == 0, z[0] == zs, v[0] == vs, 
     z'[0] == 0, v'[0] == 0}, {z, v}, {x, -1, 1}, 
    MaxStepFraction -> 1/11, 
    Method -> {"BDF", "MaxDifferenceOrder" -> 4}]
  ];

dd = Solve[{Eqn1 == 0, Eqn2 == 0}, {z''[x], v''[x]}];

Plot[Evaluate[{z''[x], v''[x]} /. dd /. sol], {x, -0.81, 0.81}]

Mathematica graphics

share|improve this answer
    
Hi :) what's the color scheme you are working with? :) –  Kuba Apr 23 at 6:19
    
@Kuba The R-Pi one. :) Or my school colors. They're the same, too. ;) –  Michael E2 Apr 23 at 10:06
    
Could you share a link? I'm not really up to date with this R-Pi turmoil :) –  Kuba Apr 23 at 10:08
    

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