| Author | Message | |||
     Svitlana Subota (photinka) New member Username: photinka Post Number: 2 Registered: 02-2007 |
I have a problem with convergence of my task. There are 2 coupled equations: epsxx=eps0*(eps_per+eps_a*(cos(Theta))^2) epsyy=eps0*(eps_per) epszz=eps0*(eps_per+eps_a*(sin(Theta))^2) epsxz=eps0*eps_a*cos(Theta)*sin(Theta) epszx=eps0*eps_a*cos(Theta)*sin(Theta) epsxy=0 epsyx=0 epsyz=0 epszy=0 D_x=(epsxx*Ex+epsxz*Ez+epsxy*Ey) D_y=(epsyy*Ey+epsyz*Ez+epsyx*Ex) D_z=(epszx*Ex+epszz*Ez+epszy*Ey) D_el=vector(D_x,D_y,D_z) EQUATIONS U: div(D_el)=0 Theta: div(grad(Theta))-eps_a*eps0/K*(Ex*cos(Theta)+Ez*sin(Theta))*(Ex*sin(Theta)-Ez*co s(Theta))=0 If I try to solve one equation, while I set another variable as constant, everything is ok. I used different values of solution method parameters, but the problem diverges. Could you advise me anything? | |||
| Robert G. Nelson (rgnelson) Moderator Username: rgnelson Post Number: 1175 Registered: 06-2003 |
Without the full script, we can only guess. One guess is that if the initial values are far from the correct solution, the newton's method tries to take jumps in Theta that cause reversals in slope. Try using SELECT CHANGELIM=0.01 to force the Newton's method to creep toward a solution. If you post or email the full script, we will be able to get a better idea of what is going wrong. | |||
| Svitlana Subota (photinka) Junior Member Username: photinka Post Number: 3 Registered: 02-2007 |
I did try it. Here is the script. SELECT VANDENBERG =on XERRLIM=0.0005 CHANGELIM =0.01 !order=3 SELECT errlim=1.e-4 COORDINATES cartesian3 VARIABLES U Theta DEFINITIONS eps0=8.854e-12 eps_per=5.3 eps_a=16.4 K=1.37E-11 U0=15 ! Applied voltage H=1.e-3 R0=1e-3 D=1*R0 epsxx=eps0*(eps_per+eps_a*(cos(Theta))^2) epsyy=eps0*(eps_per) epszz=eps0*(eps_per+eps_a*(sin(Theta))^2) epsxz=eps0*eps_a*cos(Theta)*sin(Theta) epszx=eps0*eps_a*cos(Theta)*sin(Theta) epsxy=0 epsyx=0 epsyz=0 epszy=0 D_x=(epsxx*Ex+epsxz*Ez+epsxy*Ey) D_y=(epsyy*Ey+epsyz*Ez+epsyx*Ex) D_z=(epszx*Ex+epszz*Ez+epszy*Ey) D_el=vector(D_x,D_y,D_z) INITIAL VALUES U=U0*z/H Theta=1.55*sin(pi*z/H) EQUATIONS U: div(D_el)=0 Theta: div(grad(Theta))-eps_a*eps0/K*(Ex*cos(Theta)+Ez*sin(Theta))*(Ex*sin(Theta)-Ez*co s(Theta))=0 EXTRUSION SURFACE 'Bottom' z=0 LAYER 'Hat' SURFACE 'Top' z=H BOUNDARIES Surface 'Bottom' value( U)=0 value(Theta)=0 Surface 'Top' value( U)=U0 value(Theta)=0 Region 'BOX' START(-D,-D) natural(U)=epsyz*dz(U) natural(Theta)=0 LINE TO (D,-D) natural(U)=-epsxz*dz(U) natural(Theta)=0 LINE TO (D,D) natural(U)=-epsyz*dz(U) natural(Theta)=0 LINE TO (-D,D) natural(U)=epsxz*dz(U) natural(Theta)=0 line TO CLOSE END | |||
| Robert G. Nelson (rgnelson) Moderator Username: rgnelson Post Number: 1178 Registered: 06-2003 |
1) Your script as sent has not defined Ex, Ey or Ez, so it cannot run. I presume you want Ex=dx(U), etc. 2) you have no plots or monitors. It is impossible to evaluate the solution progress without some kind of visual report. 3) Your solution for theta rides exactly at pi/2 for most of the domain. At this value, product terms of sines and cosines are ambiguous: (cos(pi/2+eps))^2 = (cos(pi/2-eps))^2, etc. I believe it is this ambiguity that prevents the Newton's method from converging. 4) This problem seems to have only z-dependence. It would be MUCH cheaper to experiment with a 1D solution until the behavior of the system is determined. 5) Given a difficult nonlinear system to solve, it is sometimes effective to use a pseudo-time relaxation system to search for the solution. This puts the evolution of the system under the control of timestep controls. 6) Given the complexity of the Jacobian derivatives of the theta equation, it may be better to hide the source term inside a SAVE, to prevent differentiation of these terms. In the attached script, I have done these things, and the solution converges quite quickly. As the attached plot shows, Theta wants to be pi/2 everywhere except where it is held to zero by boundary conditions.
 | |||
| Robert G. Nelson (rgnelson) Moderator Username: rgnelson Post Number: 1179 Registered: 06-2003 |
Here is a 3D pseudo-time script.
| |||
| Svitlana Subota (photinka) Member Username: photinka Post Number: 4 Registered: 02-2007 |
Thank you for your help. 4) This problem is only auxiliary for the one, which is unhomogeneous in 3d. 6) Isn't the result changed because the sorce contains variables theta and U? | |||
| Robert G. Nelson (rgnelson) Moderator Username: rgnelson Post Number: 1180 Registered: 06-2003 |
4) I realize that this is an abstraction preliminary to a more complex model. But as long as the problem doesn't run, it is a moot point. Once you devise a strategy that works in a 1D or 2D model, you can proceed to the real thing. It's just that 1D is hundreds of times faster than 3D, so it's more convenient for experimentation. 6) SAVE functions are re-evaluated at each timestep, so in a relaxation model like this one, the final result has incorporated fully updated values in the SAVE argument. They just don't get differentiated in forming the Jacobian matrix. | |||
| Robert G. Nelson (rgnelson) Moderator Username: rgnelson Post Number: 1181 Registered: 06-2003 |
Actually, it appears that the SAVE is not necessary. The pseudo-time relaxation seems to be stable enough without it. |
