differential equation

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

differential equation

CdeMills
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

Liam Groener-2

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

Thomas Shores-2
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);")

_______________________________________________
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

Liam Groener-2


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

Richard Crozier
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
Reply | Threaded
Open this post in threaded view
|

Re: differential equation

Joshua Stults
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

Thomas Shores-2
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