Message Boards Message Boards

0
|
6075 Views
|
6 Replies
|
2 Total Likes
View groups...
Share
Share this post:

DSolve returns Infinite expression 1/0

Posted 10 years ago
Hi all.
I'm using mathematica in order to solve some differential equations. Using DSolve the program returns me: Power::infy: "Infinite expression 1/0. encountered."
I can't really think about a possible solution about this. Here's the code:
 p = 10;
 E1 = 210000;
 R = 100;
 a = 30;
 s = 10;
 v = 0.3;
 
 t1 = p*r/2;
 t2 = p*(R^2 - r^2)/(2*r);
B = E1*s^3/(12*(1 - v^2));

sol = DSolve[{D[1/r*D[(r*w1'[r]), r], r] == t1/B,
     D[1/r*D[(r*w2'[r]), r], r] == t2/B, w1'[0] == 0, w1[a] == 0,
     w2[a] == 0, w1'[a] == w2'[a],
     w1''[a] + v*(w1'[a]/a) == w2''[a] + v*(w2'[a]/a),
     w2''[R] + v*w2'[R]/R == 0}, {w1[r], w2[r]}, r];
Can anyone help me please? Many thanks.
POSTED BY: Andrea Gazzotti
6 Replies
Posted 10 years ago
Multiplying both the left and right side of the first two equations by r^2 clears r from the denominator. Then DSolve can solve the system.
(I edited ths since first posted -- I didn't notice the independent variable in the earlier assignment to t1 and t2.)

 In[20]:= t1 = p*r/2;
 t2 = p*(R^2 - r^2)/(2*r);
 
 In[22]:= eqs2 = {r^2 D[1/r*D[(r*w1'[r]), r], r] == r^2 t1/B,
     r^2 D[1/r*D[(r*w2'[r]), r], r] == r^2 t2/B, w1'[0] == 0,
     w1[a] == 0, w2[a] == 0, w1'[a] == w2'[a],
     w1''[a] + v*(w1'[a]/a) == w2''[a] + v*(w2'[a]/a),
     w2''[R] + v*w2'[R]/R == 0} // Simplify;
 
In[23]:= DSolve[eqs2, {w1[r], w2[r]}, r]

Out[23]= {{w1[r] -> (1/(
   64 B R^2 (1 + v)))(4 a^6 p - 4 a^4 p r^2 + 3 a^4 p R^2 -
     4 a^2 p r^2 R^2 + p r^4 R^2 - 2 a^2 p R^4 + 2 p r^2 R^4 -
     4 a^6 p v + 4 a^4 p r^2 v + 11 a^4 p R^2 v - 12 a^2 p r^2 R^2 v +
      p r^4 R^2 v - 6 a^2 p R^4 v + 6 p r^2 R^4 v -
     8 a^2 p R^4 Log[a] + 8 p r^2 R^4 Log[a] - 8 a^2 p R^4 v Log[a] +
     8 p r^2 R^4 v Log[a] + 8 a^2 p R^4 Log[R] - 8 p r^2 R^4 Log[R] +
     8 a^2 p R^4 v Log[R] - 8 p r^2 R^4 v Log[R]),
  w2[r] -> (1/(
   64 B R^2 (1 + v)))(4 a^6 p - 4 a^4 p r^2 - 3 a^4 p R^2 +
     4 a^2 p r^2 R^2 - p r^4 R^2 + 6 a^2 p R^4 - 6 p r^2 R^4 -
     4 a^6 p v + 4 a^4 p r^2 v + 5 a^4 p R^2 v - 4 a^2 p r^2 R^2 v -
     p r^4 R^2 v + 2 a^2 p R^4 v - 2 p r^2 R^4 v +
     8 a^4 p R^2 Log[a] - 16 a^2 p R^4 Log[a] +
     8 a^4 p R^2 v Log[a] - 16 a^2 p R^4 v Log[a] -
     8 a^4 p R^2 Log[r] + 8 a^2 p R^4 Log[r] + 8 p r^2 R^4 Log[r] -
     8 a^4 p R^2 v Log[r] + 8 a^2 p R^4 v Log[r] +
     8 p r^2 R^4 v Log[r] + 8 a^2 p R^4 Log[R] - 8 p r^2 R^4 Log[R] +
     8 a^2 p R^4 v Log[R] - 8 p r^2 R^4 v Log[R])}}

With the values (not included above):

POSTED BY: David Keith
Since you have numerical parameters I am guessing that you are not necessarilly interested in an analytical solution.  If that is true then you can use NDSolve with a range for parameter r that avoids the potential singularity at r=0.  (I am not saying that there is indeed a singularity there, but if your range of values that you are interested in  does not include r=0, then that will be a moot question.  Also you can use NDSolve to explore whether there is indeed a singularity at r=0.)

On looking at this further I notice that there still are w/0 errors generated when using NDSolve...  stay tuned....

Ok, change your value for v to 
v = 3/10;
instead of a floating point value. Then DSolve works for you.  I am not completelysure why the use of the floating point value causes DSolve to fail, but this does create a result. 
  p = 10;
  E1 = 210000;
  R = 100;
  a = 30;
  s = 10;
  v = 0.3;
 
  t1 = p*r/2;
  t2 = p*(R^2 - r^2)/(2*r);
B = E1*s^3/(12*(1 - v^2));


sol = DSolve[{D[1/r*D[(r*w1'[r]), r], r] == t1/B,
     D[1/r*D[(r*w2'[r]), r], r] == t2/B, w1'[0] == 0, w1[a] == 0,
     w2[a] == 0, w1'[a] == w2'[a],
     w1''[a] + v*(w1'[a]/a) == w2''[a] + v*(w2'[a]/a),
     w2''[R] + v*w2'[R]/R == 0}, {w1[r], w2[r]}, r];

And this yields the following for the value of sol:
{{w1[r] -> (1/
   1600000000)(-288928800 + 309332 r^2 + 13 r^4 - 936000000 Log[30] +
     1040000 r^2 Log[30] + 936000000 Log[100] - 1040000 r^2 Log[100]),
   w2[r] -> (1/
   1600000000)(583891200 - 637068 r^2 - 13 r^4 - 1787760000 Log[30] +
     936000000 Log[100] - 1040000 r^2 Log[100] + 851760000 Log[r] +
     1040000 r^2 Log[r])}}

If you now explore your functions near r=0 you see that W2 has a singularity there:
Plot[Evaluate[{w1[r], w2[r]} /. sol[[1]]], {r, 0, 10},
PlotRange -> All]
POSTED BY: David Reiss
Posted 10 years ago
Sir you are a genius. Really thank you.
POSTED BY: Andrea Gazzotti
A genius would have figured out why the floating point value casued the problem! Perhaps this is a bug....  or perhaps Mathematica explores numerical possibilities when parameters in the differential equations are not purly symbolic/exact.  

But you are very welcome!
POSTED BY: David Reiss
Possibly rewriting the equations without r in the denominators would help.
POSTED BY: Frank Kampas
Posted 10 years ago
Any guess about how can i manage to do that?
POSTED BY: Andrea Gazzotti
Reply to this discussion
Community posts can be styled and formatted using the Markdown syntax.
Reply Preview
Attachments
Remove
or Discard

Group Abstract Group Abstract