# differential equation

 Classic List Threaded
7 messages
Reply | Threaded
Open this post in threaded view
|

## differential equation

 Hello, I'm trying to simulate an electrical circuit whose time-domain differential equation is V" + k1 V' +k2 V = a1 I' + a2 I Vs(t) - Rs I = V where V and I are the device voltage and current, respectively, and the rest are constants. Vs is the time-varying voltage source, and Rs its internal resistance. Starting from equilibrium, all values and their derivatives are set to zero. I don't see how to formulate the problem, in particular adding the source equation into the picture. LSODE: let's say x1 = V, x2 = I, x3 = V'; so we have x1' = x3 x2' = ??? => I' depends on Vs'(t) ??? x3' = a1 * ??? + a2 x2 - k1 x3 - k2 x1 How to fill the "???" DASSL: same variables res1 = x1'-x3 res2 = ... res3 = x3' - a1 x2' - a2 x2 - k1 x3 - k2 x1 In all cases, I need a relationship around I', and would like Vs(t) to be a step source, in this case Vs'(t) would become infinite ??? Any hints ? Regards Pascal
Reply | Threaded
Open this post in threaded view
|

## Re: differential equation

 On May 16, 2011, at 9:35 AM, CdeMills wrote: > Hello, > > I'm trying to simulate an electrical circuit whose time-domain differential > equation is > V" + k1 V' +k2 V = a1 I' + a2 I > Vs(t) - Rs I = V > where V and I are the device voltage and current, respectively, and the rest > are constants. Vs is the time-varying voltage source, and Rs its internal > resistance. Starting from equilibrium, all values and their derivatives are > set to zero. > > I don't see how to formulate the problem, in particular adding the source > equation into the picture. > > LSODE: let's say x1 = V, x2 = I, x3 = V'; so we have > x1' = x3 > x2' = ??? => I' depends on Vs'(t) ??? > x3' = a1 * ??? + a2 x2 - k1 x3 - k2 x1 > How to fill the "???" > > DASSL: same variables > res1 = x1'-x3 > res2 = ... > res3 = x3' - a1 x2' - a2 x2 - k1 x3 - k2 x1 > In all cases, I need a relationship around I', and would like Vs(t) to be a > step source, in this case Vs'(t) would become infinite ??? > > Any hints ? I would substitute I = (Vs-V)/Rs into your first equation and use lsode. You would then have 2 d.e.'s: x'(1) = x(2) x'(2) = aVs - k1 x(2) - b x(1) where: a = a1/Rs and b = k2+a1/Rs The initial conditions would be: x(1)=0, x(2)=0 Don't worry about the infinite derivative, the imbalance in the initial voltage should take care of that.   _______________________________________________ Help-octave mailing list [hidden email] https://mailman.cae.wisc.edu/listinfo/help-octave
Reply | Threaded
Open this post in threaded view
|

## Re: differential equation

 On 05/16/2011 04:52 PM, Liam Groener wrote: > On May 16, 2011, at 9:35 AM, CdeMills wrote: > >> Hello, >> >> I'm trying to simulate an electrical circuit whose time-domain differential >> equation is >> V" + k1 V' +k2 V = a1 I' + a2 I >> Vs(t) - Rs I = V >> where V and I are the device voltage and current, respectively, and the rest >> are constants. Vs is the time-varying voltage source, and Rs its internal >> resistance. Starting from equilibrium, all values and their derivatives are >> set to zero. >> >> I don't see how to formulate the problem, in particular adding the source >> equation into the picture. >> >> LSODE: let's say x1 = V, x2 = I, x3 = V'; so we have >> x1' = x3 >> x2' = ??? =>  I' depends on Vs'(t) ??? >> x3' = a1 * ??? + a2 x2 - k1 x3 - k2 x1 >> How to fill the "???" >> >> DASSL: same variables >> res1 = x1'-x3 >> res2 = ... >> res3 = x3' - a1 x2' - a2 x2 - k1 x3 - k2 x1 >> In all cases, I need a relationship around I', and would like Vs(t) to be a >> step source, in this case Vs'(t) would become infinite ??? >> >> Any hints ? > > I would substitute I = (Vs-V)/Rs into your first equation and use lsode. > You would then have 2 d.e.'s: > x'(1) = x(2) > x'(2) = aVs - k1 x(2) - b x(1) > where: > a = a1/Rs and b = k2+a1/Rs > > The initial conditions would be: > x(1)=0, x(2)=0 > Don't worry about the infinite derivative, the imbalance in the initial voltage should take care of that. > > _______________________________________________ > Help-octave mailing list > [hidden email] > https://mailman.cae.wisc.edu/listinfo/help-octaveInteresting.  This works.  In any case, the three variable cannot possibly work, since lsode is based on methods (ABM, BDF) that assume continuous solutions, and this is all they will generate.  However, the formula Vs = V+Rs*I guarantees that one of V or I cannot possibly be continuous at any jump discontinuity of Vs.  Thus the 3-variable approach cannot possibly work in any formulation since it would mandate that V, V' and I be continuous initially and at every interior point.   In that regard, the conditions that I and I' vanish initially make the system inconsistent if they are meant to be used as part of an initial value problem.   Your approach ignore those constraints and  banished the potential discontinuities due to the step source on the right hand side of the first order ode as simple jump discontinuities, so the problem remains well posed.  Here's a script that illustrates these points with an impulse function. Thomas Shores % DEtest.m clear global k = [1,2] # [k1,k2] global a = [2,2] # [a1,a2] global Rs = 0.5 function retval = heaviside(t) retval = (t==0)*0.5 + (t>0); end function retval = Vs(t) retval = heaviside(t-0.5)-heaviside(t-1); end function retval = Vsp(t) retval = zeros(size(t)); # ignore any dirac delta function end function retval = fcn(x,t) # here x = [V; Vp] global k; global a; global Rs; retval = [0;0]; retval(1) = x(2); retval(2) = a(1)*(Vsp(t)-x(2))/Rs + a(2)*(Vs(t)-x(1))/Rs - k(1)*x(2) - k(2)*x(1); end % main: give it a spin x0 = [0;0]; t = linspace(0,4,800); x = lsode(@fcn,x0,t); clf plot(t',x(:,1),"-1;V(t);") hold on, grid plot(t',x(:,2),"-2;V'(t);") plot(t',(Vs(t)'-x(:,1))/Rs,"-3;I(t);") plot(t,Vs(t),"-4;Vs(t);") _______________________________________________ Help-octave mailing list [hidden email] https://mailman.cae.wisc.edu/listinfo/help-octave
Reply | Threaded
Open this post in threaded view
|

## Re: differential equation

 Sent from my iPad On May 16, 2011, at 8:41 PM, Thomas Shores <[hidden email]> wrote: > On 05/16/2011 04:52 PM, Liam Groener wrote: >> On May 16, 2011, at 9:35 AM, CdeMills wrote: >> >>> Hello, >>> >>> I'm trying to simulate an electrical circuit whose time-domain differential >>> equation is >>> V" + k1 V' +k2 V = a1 I' + a2 I >>> Vs(t) - Rs I = V >>> where V and I are the device voltage and current, respectively, and the rest >>> are constants. Vs is the time-varying voltage source, and Rs its internal >>> resistance. Starting from equilibrium, all values and their derivatives are >>> set to zero. >>> >>> I don't see how to formulate the problem, in particular adding the source >>> equation into the picture. >>> >>> LSODE: let's say x1 = V, x2 = I, x3 = V'; so we have >>> x1' = x3 >>> x2' = ??? =>  I' depends on Vs'(t) ??? >>> x3' = a1 * ??? + a2 x2 - k1 x3 - k2 x1 >>> How to fill the "???" >>> >>> DASSL: same variables >>> res1 = x1'-x3 >>> res2 = ... >>> res3 = x3' - a1 x2' - a2 x2 - k1 x3 - k2 x1 >>> In all cases, I need a relationship around I', and would like Vs(t) to be a >>> step source, in this case Vs'(t) would become infinite ??? >>> >>> Any hints ? >> >> I would substitute I = (Vs-V)/Rs into your first equation and use lsode. >> You would then have 2 d.e.'s: >> x'(1) = x(2) >> x'(2) = aVs - k1 x(2) - b x(1) >> where: >> a = a1/Rs and b = k2+a1/Rs >> >> The initial conditions would be: >> x(1)=0, x(2)=0 >> Don't worry about the infinite derivative, the imbalance in the initial voltage should take care of that. >> >> _______________________________________________ >> Help-octave mailing list >> [hidden email] >> https://mailman.cae.wisc.edu/listinfo/help-octave> Interesting.  This works.  In any case, the three variable cannot possibly work, since lsode is based on methods (ABM, BDF) that assume continuous solutions, and this is all they will generate.  However, the formula Vs = V+Rs*I guarantees that one of V or I cannot possibly be continuous at any jump discontinuity of Vs.  Thus the 3-variable approach cannot possibly work in any formulation since it would mandate that V, V' and I be continuous initially and at every interior point.  In that regard, the conditions that I and I' vanish initially make the system inconsistent if they are meant to be used as part of an initial value problem.   Your approach ignore those constraints and  banished the potential discontinuities due to the step source on the right hand side of the first order ode as simple jump discontinuities, so the problem remains well posed.  Here's a script that illustrates these points with an impulse function. > > Thomas Shores > > % DEtest.m > > clear > global k = [1,2] # [k1,k2] > global a = [2,2] # [a1,a2] > global Rs = 0.5 > > function retval = heaviside(t) > retval = (t==0)*0.5 + (t>0); > end > > function retval = Vs(t) > retval = heaviside(t-0.5)-heaviside(t-1); > end > > function retval = Vsp(t) > retval = zeros(size(t)); # ignore any dirac delta function > end > > function retval = fcn(x,t) > # here x = [V; Vp] > global k; > global a; > global Rs; > > retval = [0;0]; > retval(1) = x(2); > retval(2) = a(1)*(Vsp(t)-x(2))/Rs + a(2)*(Vs(t)-x(1))/Rs - k(1)*x(2) - k(2)*x(1); > end > > % main: give it a spin > x0 = [0;0]; > t = linspace(0,4,800); > x = lsode(@fcn,x0,t); > clf > plot(t',x(:,1),"-1;V(t);") > hold on, grid > plot(t',x(:,2),"-2;V'(t);") > plot(t',(Vs(t)'-x(:,1))/Rs,"-3;I(t);") > plot(t,Vs(t),"-4;Vs(t);") > Actually, I see that I forgot the I' part of the original d.e., so back to the drawing boards. _______________________________________________ Help-octave mailing list [hidden email] https://mailman.cae.wisc.edu/listinfo/help-octave
Reply | Threaded
Open this post in threaded view
|

## Re: differential equation

 Must Vs be a perfect step source, or could you approximate it with a very steep slope? Presumably in reality it would take some finite time to go from one value to the other? You could then specify the derivative algebraically in your function something along the lines of: if        0 <= t < t1            Vs' = 0; elseif   t1 <= t < t2           Vs' = 1000; (or perhaps make Vs some quadratic function during this period?) elseif   t2 <= t < tend        Vs' = 0; Just a suggestion -- View this message in context: http://octave.1599824.n4.nabble.com/differential-equation-tp3526715p3528765.htmlSent from the Octave - General mailing list archive at Nabble.com. _______________________________________________ Help-octave mailing list [hidden email] https://mailman.cae.wisc.edu/listinfo/help-octave
Reply | Threaded
Open this post in threaded view
|

## Re: differential equation

 In reply to this post by Liam Groener-2 Oops, forgot to reply to the list. ---------- Forwarded message ---------- From: Joshua Stults <[hidden email]> Date: Tue, May 17, 2011 at 9:43 AM Subject: Re: differential equation To: Liam Groener <[hidden email]> On Tue, May 17, 2011 at 1:04 AM, Liam Groener <[hidden email]> wrote: > Actually, I see that I forgot the I' part of the original d.e., so back to the drawing boards. Extending your original approach, you can take the derivative of the Vs - Rs I = V equation and eliminate I' as well as I; then you have a single second order ode in V with a forcing term that includes Vs and Vs'.  This doesn't solve the op's problem with infinite Vs' though. -- Joshua Stults Website: variousconsequences.com Hackerspace: daytondiode.org -- Joshua Stults Website: variousconsequences.com Hackerspace: daytondiode.org _______________________________________________ Help-octave mailing list [hidden email] https://mailman.cae.wisc.edu/listinfo/help-octave
Reply | Threaded
Open this post in threaded view
|

## Re: differential equation

 In reply to this post by Richard Crozier On 05/17/2011 05:38 AM, Richard Crozier wrote: ```Must Vs be a perfect step source, or could you approximate it with a very steep slope? Presumably in reality it would take some finite time to go from one value to the other? You could then specify the derivative algebraically in your function something along the lines of: if 0 <= t < t1 Vs' = 0; elseif t1 <= t < t2 Vs' = 1000; (or perhaps make Vs some quadratic function during this period?) elseif t2 <= t < tend Vs' = 0; Just a suggestion -- View this message in context: http://octave.1599824.n4.nabble.com/differential-equation-tp3526715p3528765.html Sent from the Octave - General mailing list archive at Nabble.com. _______________________________________________ Help-octave mailing list [hidden email] https://mailman.cae.wisc.edu/listinfo/help-octave ``` I have responded to this topic twice two days ago and my posts have not yet appeared on help-octave, so let me try one more time.  This was in response to Liam Groener's response to a prior post: ```Actually, I see that I forgot the I' part of the original d.e., so back to the drawing boards. ``` Sorry if this is a repost, but my reply hasn't appeared on octave-help, so here it is (maybe again). Drats!  You're not alone.  I neglected the effect of the derivatives of Vs(t) as well, which is incorrect  So it's back to applied math 101, where the real problem with this model lies.  lsode handles jump discontinuities in the rhs just fine, and produces a continuous solution, as it must.  However, it cannot handle delta functions on the rhs, because the solution in this case *must* also suffer a discontinuity.  So the right thing to do is to restart lsode at every jump discontinuity of Vs(t) and account for the delta function by adding the jump to the new initial condition for the corresponding lhs, which in this case is V'.  So:  Your approach is correct in ignoring the initial conditions on I, and if you make the modifications I suggest it will produce a correct solution.  Where Vs(t) suffers a jump discontinuity, both I(t) and V'(t) will do the same, but V(t) remains continuous across the entire interval.  I modified my original code to yield a correct (I think) solution to the problem.  Here is the script for a sample discontinuous impulse source term Vs.  The graph is busy, so rem out the plots you don't want to see: % DEtest.m clear global k = [1,2] # [k1,k2] global a = [2,2] # [a1,a2] global Rs = 0.5 function retval = heaviside(t) retval = (t==0)*0.5 + (t>0); end function retval = Vs(t) retval = heaviside(t-0.5)-1.5*heaviside(t-1)+0.5*heaviside(t-2); end function retval = Vsp(t) retval = zeros(size(t)); % only valid off discontinuities end function retval = fcn(x,t) # here x = [V; Vp] global k; global a; global Rs; retval = [x(2); a(1)*(Vsp(t)-x(2))/Rs+a(2)*(Vs(t)-x(1))/Rs-k(1)*x(2)-k(2)*x(1)]; end % main: give it a spin bdy = [0.5,1; 1,-1.5; 2,0.5; 4,0]; # [point, jump] pairs for lsode up to right bdy xsoln = x0 = [0,0]; tsoln = 0; for j = 1:rows(bdy)   t = linspace(tsoln(end),bdy(j,1),100);   x = lsode(@fcn,x0,t);   tsoln = [tsoln,t(2:end)];   xsoln = [xsoln;x(2:end,:)];   x0 = (xsoln(end,:)+[0,bdy(j,2)])'; end clf plot(tsoln',xsoln(:,1),"-1;V(t);") hold on, grid plot(tsoln',xsoln(:,2),"-2;V'(t);") plot(tsoln',(Vs(tsoln)'-xsoln(:,1))/Rs,"-3;I(t);") plot(tsoln,Vs(tsoln),"-4;Vs(t);") _______________________________________________ Help-octave mailing list [hidden email] https://mailman.cae.wisc.edu/listinfo/help-octave