@dharr 

Thank you, I have since tried some direct dsolve solutions with the same results as convertsys.And then I'm going to go straight to dsolve.

@C_R 

Sorry, I saw this error after not able to solve, originally thought is not called the Jacobi function package, but after the call still appeared this error

@C_R 

Thanks for your guidance and understanding that global variables cannot be directly replaced with local variables. I saw above that you said there was a problem with this Jacobian matrix. Can you elaborate?

@C_R 

@C_R 

Thanks, (tan(tehta))^2=-1 has no Angle correspondence. Combined with the actual situation of the problem, it is found that theta has little influence on the dynamic response. After solving a and c, I can estimate the result。Sincerely thank you for your approach, I learned a new method to solve this question.I find my maple ability is too poor, I would like to ask you what maple books recommend?

@C_R 

Learn the code you wrote and feel pretty good, but find that plugging a1 and c1 into secular7[1] is not the same as solving theta3 in secular7[3].

@Carl Love 

Thank you for your reply, using your method can indeed jump the singularity and continue to solve, but there are still many singularity points later, and the distance between the two singular points is getting closer and closer, I would like to ask why?

solve826.mw

@mmcdara 

Thank you for reminding me. I did not realize that gamma has a built-in function before. I have modified gamma

solu := dsolve(eval({ics, couple[]}, gamma = gamma__1), [eta(t), varphi(t), V(t)], output = procedurelist, numeric, maxfun = -1)

solve826.mw

@Carl Love When this method is used, the same command is used, but the error will always be reported. Why is this?

solve826.mw

@C_R I think it is more reasonable to solve the problem in this way

@C_R  Is this equation unreasonable?

@C_R Is that how it's handled?

@mmcdara solve821.mw

solve821.mw
 

Download solve821.mw
 

Download solve821.mwsolve821.mw

@Rouben Rostamian  

restart;
secular7 := -3.461584716*10^11*a__1*((1.024337843*a__1^4 + (4.073013529*c__1^2 + 29.75257330)*a__1^2 + c__1^4 + 25.6247057300000*c__1^2 - 9.31744223200000*10^(-9))*cos(theta__3) + (0.1403469826*a__1^2 + 0.1233422711*c__1^2 - 2.396674825)*sin(theta__3)), -3.752578062*10^10*((1.000000000*c__1^4 + (12.21904058*a__1^2 + 13.28028686)*c__1^2 - 110.1848865*a__1^2 + 9.219040578*a__1^4 + 1.057632194*10^(-8))*cos(theta__4) + (-0.1226429622*c__1^2 - 0.5582050515*a__1^2 + 4.766172985)*sin(theta__4))*c__1, -3.461584716*10^11*a__1*((1.024337843*a__1^4 + (4.073013529*c__1^2 + 29.75257330)*a__1^2 + c__1^4 + 25.6247057300000*c__1^2 - 9.31744223200000*10^(-9))*sin(theta__3) + (-0.1403469826*a__1^2 - 0.1233422711*c__1^2 + 2.396674825)*cos(theta__3)), -3.752578062*10^10*((1.000000000*c__1^4 + (12.21904058*a__1^2 + 13.28028686)*c__1^2 - 110.1848865*a__1^2 + 9.219040578*a__1^4 + 1.057632194*10^(-8))*sin(theta__4) + (0.1226429622*c__1^2 + 0.5582050515*a__1^2 - 4.766172985)*cos(theta__4))*c__1;
                           11      //                4
secular7 := -3.461584716 10   a__1 \\1.024337843 a__1 

     /                2              \     2       4
   + \4.073013529 c__1  + 29.75257330/ a__1  + c__1 

                          2                 -9\                 
   + 25.6247057300000 c__1  - 9.317442232 10  / cos(theta__3) + 

  /                 2                    2              \     
  \0.1403469826 a__1  + 0.1233422711 c__1  - 2.396674825/ sin(

           \                 10 //                4
  theta__3)/, -3.752578062 10   \\1.000000000 c__1 

     /                2              \     2                   2
   + \12.21904058 a__1  + 13.28028686/ c__1  - 110.1848865 a__1 

                     4                 -8\                 /
   + 9.219040578 a__1  + 1.057632194 10  / cos(theta__4) + \
                  2                    2              \          
-0.1226429622 c__1  - 0.5582050515 a__1  + 4.766172985/ sin(thet\

       \                      11      //                4
  a__4)/ c__1, -3.461584716 10   a__1 \\1.024337843 a__1 

     /                2              \     2       4
   + \4.073013529 c__1  + 29.75257330/ a__1  + c__1 

                          2                 -9\                 /
   + 25.6247057300000 c__1  - 9.317442232 10  / sin(theta__3) + \
                  2                    2              \          
-0.1403469826 a__1  - 0.1233422711 c__1  + 2.396674825/ cos(thet\

       \                 10 //                4
  a__3)/, -3.752578062 10   \\1.000000000 c__1 

     /                2              \     2                   2
   + \12.21904058 a__1  + 13.28028686/ c__1  - 110.1848865 a__1 

                     4                 -8\                 
   + 9.219040578 a__1  + 1.057632194 10  / sin(theta__4) + 

  /                 2                    2              \     
  \0.1226429622 c__1  + 0.5582050515 a__1  - 4.766172985/ cos(

           \     
  theta__4)/ c__1

with(RealDomain);
solv := fsolve({secular7}, {a__1, c__1, theta__3, theta__4});
    solv := {a__1 = 0., c__1 = 0., theta__3 = -8.133518273, 

      theta__4 = 6.228573980}

@janhardo I'm sorry, I didn't understand what this comment meant? Could you elaborate?