---------- JoTz wrote: ----------
As I'm having a slight problem with the email to this group , I I shall attempt to send to you. I interpret this to mean that I should post this to the [Fractint] list. Done.
I have some concern about the email systems and / or list serve software involved reformatting the layout of the text... I don't think I've seen characters thrown away -- just lines re-formatted so they become less easy to read. Re what JoTz's equations are modeling:
The circuit [simulation] supposedly generates a rectified sine wave [over time]; . . . The equation is a second order , very stiff non linear differential equation...
- Hal Lane ######################## # hallane@earthlink.net ######################## From: sciwise@ihug.co.nz [mailto:sciwise@ihug.co.nz] Sent: Sunday, January 10, 2016 10:15 PM To: Hal Lane <hallane@earthlink.net> Subject: RE: [Fractint] Tedious formula Hi Hal, As I'm having a slight problem with the email to this group , I I shall attempt to send to you. I might be able to replace the exp with temp*temp. To be really accurate exp() is what's used , however to get a general idea as to what's going on your recommendation has merit. I tried setting the initial (x,y) to (0,0) , I managed to get a few more iterations , however the overflow problems remain. Setting (p,h) to a small value provides a few more iterations, however the output isn't quite what I expected ; I may need to use a circuit simulator , like TopSpice to determine what the output is meant to look like. The equation is a second order , very stiff non linear differential equation , I haven't been able to find an analytical solution to this upon the internet. You're right about the reuse of x and y , however the first (x,y) values are outside the main loop also this type of coding isn't different from the Volterra Lotka formula file that comes standard with Fractint. I'm fairly certain that alteration of (x,y) within the main loop is correct. Nevertheless I shall try your recommendation. The q value is reused , within the main loop , however I'm not expecting any difficulty in doing so ; however I can also change that. The variable t1 is an independent time variable , used when calculating the value of sin() ; the input function. I'm now pondering how I might do a two loop iteration within Fractint ; can't figure out how it might be done. The circuit supposedly generates a rectified sine wave ; however I really need to look at the circuit simulator results. This is all a bit of a jigsaw puzzle. I'm including my formula file for the diode model as this illustrates how I dealt with the exp() term in a simpler instance. ======== START PARAMETER FILE======================================= Diode { reset=1821 type=formula formulafile=c:\winf1821\alndgh.frm formulaname=LnDioGH corners=-0.16588419/0.25978091/-0.15031315/0.1691023 float=y maxiter=512 logmap=yes colors=00000j0j0jj0j00j0j0jjkkkkrkdmwW000W0WW000WW0W0WWFNF<5>zNFFVF<5>zV\ FFbF<5>zbFFjF<5>zjFFrF<5>zrFFzF<5>zzFFFN<5>zFNFNN<5>zNNFVN<5>zVNFbN<5>zb\ NFjN<5>zjNFrN<5>zrNFzN<5>zzNFFV<5>zFVFNV<5>zNVFVV<5>zVVFbV<5>zbVFjV<5>zj\ VFrV<5>zrVFzV<5>zzVFFb<5>zFbFNb<5>zNbFVb<5>zVbFbb<5>zbbFjb<5>zjbFrb<5>zr\ bFzb<5>zzbFFj<5>zFjFNj<5>zNjFVj<5>zVjFbj<5>zbjFjj<5>zjjFrj<4>rrjzywccdWW\ Wz000z0zz000zz0z0zzzzz } comment { From : sciwise@ihug.co.nz July 2009 copyright c 2009 This is a model of a diode , as an electronic component . It utilizes Newton's method to solve a non-linear equation ; known as the ideal diode or Shockley equation . The very basic form of this equation is : i = Ids*(exp(Vs/(n*Vt))-1) where , Q = 1.6*(10^-19) Electronic Charge . K = 1.38*(10^-23) Boltzman's constant . T = Absolute Temperature , in Kelvins . T = 273 [K] + degreesC For degreesC = 25 , Vt ~= 0.025 [v] Ids = Saturation current . Vs = Applied , source , voltage . n = Ideality factor = ( 1 to 2 ) . Vt = K * T / Q Thermionic Voltage . A more realistic model of a diode will include the finite ohmic resistance of the leads . A parallel conductance is also modeled , as this speeds the convergence of Newton's method . Therefore the circuit we're modelling looks like this . i --> ----------------------------, ' | | -------------------, | | .---. Vs - D --- | | Gm --- /_\ |___| | |__________________| | | | .---. | | | Rs | |___| '----------------------------' Where D is represented by the ideal diode equation . Now the total current , i , flowing through this circuit is : i = Ids*(exp((Vs-i*Rs)/(n*Vt))-1) + Gm*(Vs-i*Rs) This equation needs to be solved in terms of i and be put in a form that can be readily and accurately solved through Newton's method . The steps involved are : 1) Let log(x) = (Vs-i*Rs)/(n*Vt) 2) Solve log(x) = (Vs-i*Rs)/(n*Vt) in terms of i . i = -(n*Vt*log(x)-Vs)/Rs 3) Substitute the expression for i , from (2) into i = Ids*(exp((Vs-i*Rs)/(n*Vt))-1) + Gm*(Vs-i*Rs) -(n*Vt*log(x)-Vs)/Rs = Ids*(x-1) + Gm*(Vs--((n*Vt*log(x)-Vs)/Rs)*Rs) -(n*Vt*log(x)-Vs)/Rs = Ids*(x-1) + Gm*((n*Vt*log(x)/Rs)*Rs) (Vs-n*Vt*log(x))/Rs = Ids*(x-1) + Gm*(n*Vt*log(x)) 4) Solve (3) for x , x = -((Gm*n*Rs*Vt+n*Vt)*log(x)-Vs-Is*Rs)/(Is*Rs),Is*Rs,1); 5) Find the derivative of (4) in terms of x . 1 = -(Gm*n*Rs*Vt+n*Vt)/(Is*Rs*x) 6) Now Newton's iteration formula is : x = x - f(x) / f '(x) , where f '(x) == d f(x) / dx . 7) For this example : f(x) = x = -((Gm*n*Rs*Vt+n*Vt)*log(x)-Vs-Is*Rs)/(Is*Rs),Is*Rs,1) and f '(x) = 1 = -(Gm*n*Rs*Vt+n*Vt)/(Is*Rs*x) 8) Therefore Newton's iteration formula for this example is : x = x - ( -((Gm*n*Rs*Vt+n*Vt)*log(x)-Vs-Is*Rs)/(Is*Rs),Is*Rs,1) - x )/( -(Gm*n*Rs*Vt+n*Vt)/(Is*Rs*x)-1) 9) Simplifying (8) x = x -( (-Ids*Rs-Vs)*x + Ids*Rs*x*x + (n+Gm*n*Rs)*Vt*x*log(x) ) / ( (n + Gm*n*Rs)*Vt+Ids*Rs*x ) Substitute u = Gm*n*Rs*Vt+n*Vt vr = Ids*Rs x = x-( -(Vs+vr)*x + vr*x*x + u*x*log(x) )/(vr*x+u)) In this logarithmic format overflow , and possibly underflow , is avoided . 10 ) We need to be able to convert x back to i , to do this we use the inverse of the original substitution : log(x) = (Vs-i*Rs)/(n*Vt) Solving for i , i = - ( log(x)*(n*Vt) - Vs ) / Rs Substitute Pt = n*Vt i = - ( log(x)*Pt - Vs ) / Rs Between successive iterations of Newton's formula the previous value of i is compared with the present, if the absolute difference is less than a predefined tolerance , Pd , then the iteration stops . Vs is a complex voltage source , possibly equivalent to a multi phase voltage source . This is a basic diode model and doesn't include the effects of diode junction capacitance . For small voltages a purely exponential representation of the diode model may converge quicker than the logarithmic and still remain accurate . The logarithmic model has been implemented as a BASIC language program and compared with data from a SPICE simulation , for a real voltage source , there is good agreement . For a Complex voltage source , the few samples generated from the BASIC program are in good agreement with those from MAXIMA CAS . } frm:LnDioGH{ ; ; ; Complex Voltage - Complex Current ; Diode Characteristics . ; ; for 1N4148 silicon type . ; Level 1 model . ; ; at 25 degrees C ; ; Pd = (10^(-13)), Gm = 9.4*(10^(-13)), Ids = 2.52*(10^(-09)), Rs = 0.968, n = 1.752, Vt = 0.025, u = Gm*n*Rs*Vt+n*Vt, Pt = n*Vt, vr = Ids*Rs, Vs = Pixel, x = (0.1,0), x = x-((u*x*log(x)+vr*x*x-(Vs+vr)*x)/(vr*x+u)), y = -(Pt*log(x)-Vs)/Rs, tol = (1,0): u = y; x = x-((u*x*log(x)+vr*x*x-(Vs+vr)*x)/(vr*x+u)); y = -(Pt*log(x)-Vs)/Rs, tol = |u-y|, tol > Pd } comment { The floating point option needs to be selected , 512 iterations , or more , should be selected . Also a logarithmic palette . Minibrot structures aren't immediately apparent at the default parameter settings . A parameter setting of : xmin - 0.005770310598 xmax - 0.005666943901 ymin 0.000516417395 ymax 0.000594104392 shows embedded minibrots and spirals , weather this is relevant in a real diode is unknown . You might be able to locate minibrots elsewhere on the map . } ========END PARAMETER FILE================== ---------------------------------------------------------------------- On 10/01/2016 17:08, Hal Lane wrote: ----------------------------------------------------------------------
Your tidy up is good... Thanks. ...I'm using maxima cas... Very good!
Here are a few more wild and crazy debugging ideas: [exp] grows quite large after one or two iterations... Here's an odd idea. Could you temporarily replace the exp with a function that doesn't grow so quickly -- just to get some results to look at? Perhaps a power? So, instead of: q = exp( s*a1*x ) you could try something like: temp = s*a1*x q = temp * temp or: temp = s*a1*x q = temp * temp * temp * temp q would then grow (or shrink for "temp" values less than one,) but not in such a violent way. Using even numbers of "temp" multiplied together always makes its sign always be positive -- like the result of an exp(). What does the physics of the RLC circuit say about the use of the exp() ? I'm also wondering if the initial x and y are meant to be 0 and 0 and v = (x,y) . That's certainly easy to try -- but the significance of doing this ought to be thought about, as well. Here's an odd observation -- x and y get used for two different purposes: z = pixel x = real(z) y = imag(z) and also: c = x+h*(u+b) d = y + h*(w + s*a1*b*b + s*b2*q*b + s*b3*b + b2*q*v + b1*q + s*b1) x = c y = d z = x*x + y*y Might the 2nd place where x and y are used need to be xtemp and ytemp? It looks like the values of x and y possibly needing to be constant during the calculation of a pixel might be getting compromised by their reuse -- unless that re-usage is intentional, Also, I note that q is used in two different calculations: This re-usage may or may not be intentional... q = exp(s*a1*x ) . . . q = exp(s*a1*a ) Note that the color map could be useful for debugging when you get to the stage where you're looking at resultant images. Perhaps you could have a map that has big color changes in adjacent lower values to be able to better see what is going on in the initial iterations. Is the image supposed to somehow visually depict the ringing of a diode-clamped RLC circuit? - Hal Lane ######################## # hallane@earthlink.net ######################## ------------------------Original Message------------------------ From: Fractint [mailto:fractint-bounces@mailman.xmission.com] On Behalf Of sciwise@ihug.co.nz Sent: Saturday, January 9, 2016 5:50 PM To: Fractint and General Fractals Discussion <fractint@mailman.xmission.com> Subject: Re: [Fractint] Tedious formula ------------------------------------------------------------------------- Your tidy up is good , for some reason I was having difficulty with negative exponents. Presently I'm using maxima cas to look at the values that exp(q) takes , and indeed it grows quite large after one or two iterations; to such an extent that maxima runs out of digits when printing. Without going into the details of the underlying equations , I encountered a similar difficulty a number of years ago when I wrote LnDioGH.frm ; do a web search you should be able to find the code. I'm also wondering if the initial x and y are meant to be 0 and 0 and v = (x,y) . ============================================= On 10/01/2016 03:52, Jonathan Osuch wrote: On Sat, 2016-01-09 at 16:28 +1300, sciwise@ihug.co.nz [1]wrote: ============================================= Difficulties remain , however this is a slightly improved version. For some reason the loop is exited after the first iteration ; I'm starting to think that the exponential function exp(q) might be the cause. It depends upon what you are trying to iterate during each loop. Take a look at the following. I cleaned up the code a little, but it still doesn't work. DiodeRLC float=y { f = 1/real(p1) z = pixel x = real(z) y = imag(z) p = 0.39757 h = 0.71743 s = -1 R = 1 L =0.001 C1 = 10^(-6) KT = ((1.3806485279) * 10^(-23)) * 270 twopi = 2 * pi q1 = (1.602176620898) * 10^(-19) Id = 10^(-12) a1 = q1/(KT) b1 = (KT)/(q1*C1*L) b2 = (KT)/(Id*q1*L) b3 = R/L t1 = 0 : v = 100*sin(twopi*t1) q = exp( s*a1*x ) u = y w = s*a1*y*y + s*b2*q*y + s*b3*y + b2*q*v + b1*q + s*b1 a = x+p*u b = y+p*w q = exp( s*a1*a ) c = x+h*(u+b) d = y + h*(w + s*a1*b*b + s*b2*q*b + s*b3*b + b2*q*v + b1*q + s*b1) x = c y = d t1 = t1 + f z = x*x + y*y t1 < maxit ; loops until t1 is >= maxit } Jonathan Links: ------ [1] mailto:sciwise@ihug.co.nz [2] mailto:Fractint@mailman.xmission.com [3] https://mailman.xmission.com/cgi-bin/mailman/listinfo/fractint _______________________________________________ Fractint mailing list Fractint@mailman.xmission.com https://mailman.xmission.com/cgi-bin/mailman/listinfo/fractint --- This email has been checked for viruses by Avast antivirus software. https://www.avast.com/antivirus --- --- This email has been checked for viruses by Avast antivirus software. https://www.avast.com/antivirus
Sciwise, My apologies for addressing you with JoTz nickname. JoTz, My apologies for making such free use of your nickname. Blame it on the... what was I saying, I forget? - Hal Lane ######################## # hallane@earthlink.net ######################## --- This email has been checked for viruses by Avast antivirus software. https://www.avast.com/antivirus
participants (1)
-
Hal Lane