0
$\begingroup$

I am trying to solve the following linear system consisting of four unknowns.

a = 1
b = 1.1
c = 1
w = 100
p = 4500
phi[z_, r_] := 
 a5 (8 z^5 - 40 r^2 z^3 + 15 r^4 z) + b5 (2 z^5 - r^2 z^3 - 3 r^4 z) +
   a1 (2 z^4 + z^2 r^2 - 
     r^4) Log[(Sqrt[z^2 + r^2] + z)/(Sqrt[z^2 + r^2] - z)] - 
  a1 6 z/Sqrt[z^2 + r^2] + 
  b1 (r^4 + 2 r^2 z^2 + 
     z^4) Log[(Sqrt[z^2 + r^2] + z)/(Sqrt[z^2 + r^2] - z)]
sr[z_, r_] = 
 Simplify[D[
    v (D[phi[z, r], {r, 2}] + D[phi[z, r], r]/r + 
        D[phi[z, r], {z, 2}]) - D[phi[z, r], {r, 2}], z] - 
   p w^2 r^2/3 - p w^2 (1 + 2 v) (1 + v) z^2/(6 v (1 - v))]
sz[z_, r_] = 
 Simplify[D[(2 - v) (D[phi[z, r], {r, 2}] + D[phi[z, r], r]/r + 
        D[phi[z, r], {z, 2}]) - D[phi[z, r], {z, 2}], z] + 
   p w^2 (1 + 3 v) r^2/(6 v)]
trz[z_, r_] = 
 Simplify[D[(1 - v) (D[phi[z, r], {r, 2}] + D[phi[z, r], r]/r + 
       D[phi[z, r], {z, 2}]) - D[phi[z, r], {z, 2}], r]]
NSolve[Integrate[sr[z, a], {z, -c, c}] == 0 && 
  Integrate[sr[z, b], {z, -c, c}] == 0 && 
  Integrate[trz[c, r] r, {r, a, b}] == 0 && 
  Integrate[sz[c, r] r, {r, a, b}] == 0, {a1, b1, a5, b5}]

I remember last week I managed to get a set of numerical answers, but this week Mathematica reports "Power::infy: Infinite expression 1/0. encountered." and keeps running without giving me a result. I don't think there is any infinity point in the function. Can someone tell me how to revise the code to get the answer? Thank you!

Improve this question
$\endgroup$
  • $\begingroup$ I forget to add that I have defined v to be 0.3. $\endgroup$ – frogjerry Jun 26 '19 at 2:35
0
$\begingroup$

I made a slight modification to your code and it gives an answer on 12.0.0 + OSX 10.14.5. You may want to make sure that you have nothing defined for any of the internal variables using ClearAll[]. 

ClearAll[z, v, n, a1, b1, a5, b5];
a = 1;
b = 1.1;
c = 1;
w = 100;
p = 4500;
phi[z_, r_] := 
  a5 (8 z^5 - 40 r^2 z^3 + 15 r^4 z) + 
  b5 (2 z^5 - r^2 z^3 - 3 r^4 z) + 
    a1 (2 z^4 + z^2 r^2 - 
      r^4) Log[(Sqrt[z^2 + r^2] + z)/(Sqrt[z^2 + r^2] - z)] - 
   a1 6 z/Sqrt[z^2 + r^2] + 
   b1 (r^4 + 2 r^2 z^2 + 
      z^4) Log[(Sqrt[z^2 + r^2] + z)/(Sqrt[z^2 + r^2] - z)];
sr[z_, r_] = 
  Simplify[D[
     v (D[phi[z, r], {r, 2}] + D[phi[z, r], r]/r + 
         D[phi[z, r], {z, 2}]) - D[phi[z, r], {r, 2}], z] - 
     p w^2 r^2/3 - p w^2 (1 + 2 v) (1 + v) z^2/(6 v (1 - v))];
sz[z_, r_] = 
  Simplify[D[(2 - v) (D[phi[z, r], {r, 2}] + D[phi[z, r], r]/r + 
         D[phi[z, r], {z, 2}]) - D[phi[z, r], {z, 2}], z] + 
    p w^2 (1 + 3 v) r^2/(6 v)];
trz[z_, r_] = 
  Simplify[D[(1 - v) (D[phi[z, r], {r, 2}] + D[phi[z, r], r]/r + 
        D[phi[z, r], {z, 2}]) - D[phi[z, r], {z, 2}], r]];
NSolve[Integrate[sr[z, a], {z, -c, c}] == 0 && 
  Integrate[sr[z, b], {z, -c, c}] == 0 && 
  Integrate[trz[c, r] r, {r, a, b}] == 0 && 
  Integrate[sz[c, r] r, {r, a, b}] == 0, {a1, b1, a5, b5}]

(* {{a1 -> -((1. (((-267024. + 
           1.95221*10^-10 I) + (21216.8 - 
            9.48377*10^-10 I) v + (384367. + 
            1.45108*10^-9 I) v^2 - (79069.4 + 
            7.10136*10^-10 I) v^3) (((331.813 + 
              8.65498*10^-13 I) - (3.41061*10^-13 + 
               3.38673*10^-13 I) v) (0. + (-10.592 + 
               11.916 v) (5.13*10^7/(-1. + v) + 
               5.*10^6/((-1. + v) v) - (
               2.63*10^7 v)/(-1. + v))) - 
         1. (0. + ((2.61056*10^6 - 1.16415*10^-10 I) + (
               870187. - 2.91038*10^-11 I)/v) (-10.592 + 
               11.916 v)) ((-10.592 + 
               11.916 v) (275.6/(-1. + v) - (
               275.6 v)/(-1. + v)) - 
            52.96 (-(91.12/(-1. + v)) + (149.8 v)/(-1. + v) - (
               58.68 v^2)/(-1. + v)))) - 
      1. (-1. (-1. ((-1.75455 - 
                1.06581*10^-14 I) - (5.07465 + 
                7.88861*10^-31 I) v) (1.37707 - 
               1.67258 v) + ((16.6258 + 
                0. I) - (13.2092 + 
                2.84217*10^-14 I) v) (-10.592 + 
               11.916 v)) ((-10.592 + 
               11.916 v) (275.6/(-1. + v) - (
               275.6 v)/(-1. + v)) - 
            52.96 (-(91.12/(-1. + v)) + (149.8 v)/(-1. + v) - (
               58.68 v^2)/(-1. + v))) + ((331.813 + 
              8.65498*10^-13 I) - (3.41061*10^-13 + 
               3.38673*10^-13 I) v) (-1. (1.37707 - 
               1.67258 v) (-(91.12/(-1. + v)) + (
               149.8 v)/(-1. + v) - (
               58.68 v^2)/(-1. + v)) + (-10.592 + 
               11.916 v) (29.6561/(-1. + v) - (
               173.856 v)/(-1. + v) + (
               144.2 v^2)/(-1. + v)))) (-1. (-52.96 (76. - 
               36. v) - 
            200. (-10.592 + 
               11.916 v)) (0. + ((2.61056*10^6 - 
                1.16415*10^-10 I) + (870187. - 2.91038*10^-11 I)/
               v) (-10.592 + 11.916 v)) + ((331.813 + 
              8.65498*10^-13 I) - (3.41061*10^-13 + 
               3.38673*10^-13 I) v) (0. + (-10.592 + 
               11.916 v) (-3.*10^7 + (
               5.*10^6 (1. + v) (1. + 2. v))/((-1. + v) v))))))/((
    2.15005*10^8 + 8.27372*10^-7 I)/(-1. + 
     v) - ((2.30422*10^9 + 7.2686*10^-6 I) v)/(-1. + 
     v) + ((9.07026*10^9 + 0.0000363434 I) v^2)/(-1. + 
     v) - ((1.5815*10^10 + 0.0000347323 I) v^3)/(-1. + 
     v) + ((1.23808*10^10 - 0.0000679472 I) v^4)/(-1. + 
     v) - ((3.3374*10^9 - 0.000127206 I) v^5)/(-1. + 
     v) - ((2.0946*10^8 + 0.0000546427 I) v^6)/(-1. + 
     v) + ((0. + 2.13791*10^-7 I) v^7)/(-1. + v))), 
b1 -> -(((4.23724*10^6 - 
      3.97001*10^6 I) ((-2.68178 + 
        2.65827 I) + (1.52727*10^16 - 
         1.63007*10^16 I) v - (3.42115*10^16 - 
         3.65144*10^16 I) v^2 - (5.0757*10^16 - 
         5.41736*10^16 I) v^3 + (1.95778*10^17 - 
         2.08956*10^17 I) v^4 - (1.10133*10^17 - 
         1.17546*10^17 I) v^5 - (1.58126*10^17 - 
         1.6877*10^17 I) v^6 + (2.30055*10^17 - 
         2.4554*10^17 I) v^7 - (1.02746*10^17 - 
         1.09662*10^17 I) v^8 + (1.5357*10^16 - 
         1.63907*10^16 I) v^9 - (4.88834*10^14 - 
         5.21739*10^14 I) v^10 + (1. + 0. I) v^11))/((-1. + 
      v) v ((3.37708 - 
        3.27992*10^-14 I) - (0.268331 - 
         1.44042*10^-14 I) v - (4.86114 - 2.53068*10^-14 I) v^2 + 
      1. v^3) ((3.87001 - 
        1.00568*10^15 I) - (33.9987 - 
         1.07779*10^16 I) v + (169.995 - 
         4.24259*10^16 I) v^2 - (162.46 - 
         7.39741*10^16 I) v^3 - (317.821 + 
         5.79108*10^16 I) v^4 + (595.004 + 
         1.56106*10^16 I) v^5 - (255.59 - 
         9.79744*10^14 I) v^6 + (1. + 0. I) v^7))), 
a5 -> -(((1.96393*10^-7 - 
      2.24393*10^-8 I) ((1.7977*10^28 - 

        1.57338*10^29 I) + (2.1689*10^30 - 
         1.89826*10^31 I) v - (1.10834*10^31 - 
         9.70043*10^31 I) v^2 + (2.12511*10^31 - 
         1.85993*10^32 I) v^3 - (1.52482*10^31 - 
         1.33455*10^32 I) v^4 - (1.08663*10^31 - 
         9.51038*10^31 I) v^5 + (3.66366*10^31 - 
         3.20651*10^32 I) v^6 - (4.01472*10^31 - 
         3.51376*10^32 I) v^7 + (2.31305*10^31 - 
         2.02443*10^32 I) v^8 - (6.28849*10^30 - 
         5.50381*10^31 I) v^9 + (3.8131*10^29 - 
         3.33729*10^30 I) v^10 + (4.71739*10^28 - 
         4.12875*10^29 I) v^11 - (9.53322*10^14 - 
         8.7082*10^14 I) v^12 + (1. + 0. I) v^13))/(((331.813 + 
        8.65498*10^-13 I) - (3.41061*10^-13 + 
         3.38673*10^-13 I) v) (-1. + 
      v) v ((3.37708 - 
        3.27992*10^-14 I) - (0.268331 - 
         1.44042*10^-14 I) v - (4.86114 - 2.53068*10^-14 I) v^2 + 
      1. v^3) ((3.87001 - 
        1.00568*10^15 I) - (33.9987 - 
         1.07779*10^16 I) v + (169.995 - 
         4.24259*10^16 I) v^2 - (162.46 - 
         7.39741*10^16 I) v^3 - (317.821 + 
         5.79108*10^16 I) v^4 + (595.004 + 
         1.56106*10^16 I) v^5 - (255.59 - 
         9.79744*10^14 I) v^6 + (1. + 0. I) v^7))), 
b5 -> ((2.58987*10^-6 - 
    2.11676*10^-8 I) ((4.39031*10^26 - 
      5.37159*10^28 I) + (9.38932*10^28 - 
       1.14879*10^31 I) v - (4.26607*10^29 - 
       5.21957*10^31 I) v^2 + (4.98212*10^29 - 
       6.09567*10^31 I) v^3 + (4.92008*10^29 - 
       6.01977*10^31 I) v^4 - (1.55503*10^30 - 
       1.9026*10^32 I) v^5 + (9.13331*10^29 - 
       1.11747*10^32 I) v^6 + (5.72367*10^29 - 
       7.00296*10^31 I) v^7 - (9.48111*10^29 - 
       1.16002*10^32 I) v^8 + (4.28735*10^29 - 
       5.24561*10^31 I) v^9 - (7.30931*10^28 - 
       8.94301*10^30 I) v^10 + (3.86213*10^27 - 
       4.72535*10^29 I) v^11 - (9.72616*10^14 - 
       9.71795*10^14 I) v^12 + (1. + 0. I) v^13))/(((331.813 + 
      8.65498*10^-13 I) - (3.41061*10^-13 + 
       3.38673*10^-13 I) v) (-1. + v) v (-0.888889 + 
    1. v) ((3.37708 - 
      3.27992*10^-14 I) - (0.268331 - 
       1.44042*10^-14 I) v - (4.86114 - 2.53068*10^-14 I) v^2 + 
    1. v^3) ((3.87001 - 
      1.00568*10^15 I) - (33.9987 - 
       1.07779*10^16 I) v + (169.995 - 
       4.24259*10^16 I) v^2 - (162.46 - 
       7.39741*10^16 I) v^3 - (317.821 + 
       5.79108*10^16 I) v^4 + (595.004 + 
       1.56106*10^16 I) v^5 - (255.59 - 
       9.79744*10^14 I) v^6 + (1. + 0. I) v^7))}} *)
Improve this answer
$\endgroup$
  • $\begingroup$ I forget to add one line of code where I define v to be 0.3. I have also noticed that while Mathematica is able to give me an answer when v is not defined, it reports "Power::infy: Infinite expression 1/0. encountered" when v is defined. Why is this happening? $\endgroup$ – frogjerry Jun 26 '19 at 2:45
  • $\begingroup$ Did you click on the stack trace? I ask because I get the same divide by zero problem and the place where it is failing is this: [Integral][Integral]Log[( z + Sqrt[1 + z^2])/(-z + Sqrt[1 + z^2])]/((0. + 1. I) - 1. z)^4 [DifferentialD]z [DifferentialD]z The stack trace may be found by clicking on the button with the 3 dots just to the left of the reported error message. And unfortunately, I don't know why it is failing. $\endgroup$ – Mark R Jun 26 '19 at 5:39
  • $\begingroup$ Problems disappear on my MMA version 8.0, when using infinite precision parameters v = 3/10, b = 11/10 avoiding precision problems. (By the way, use SetDelayed := only when necessary, here with phi[z_,r_]:= ). $\endgroup$ – Akku14 Mar 21 '20 at 21:15

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.