MATLAB ODE45 question

I have the following differential equation:

[eqn] LQ'' + RQ' + \frac{Q}{C} = 230 \sqrt{2} \sin(\omega \pi) [/eqn]

I want to know if I plug this into a function correctly. I put [math]Q' = I[/math]. Then we get:

[eqn] LI' + RI + \frac{Q}{C} = 230 \sqrt{2} \sin( \omega \pi ) [/eqn]

If we then define the vector [math] \mathbf{y} [/math] as [math]y_1 = I[/math] and [math]y_2 = Q[/math]. We get:

[eqn] y_1' = \frac{230 \sqrt{2} \sin ( \omega \pi )}{L} - \frac{R y_1}{L} - \frac{ y_2 }{LC} [/eqn]

[eqn] y_2' = y_1 [/eqn]

Plugging this into a function file we would get:

cont'd in next post.

This is my function file:

--------------------------------------------------------

function dydt = ode5( t,y,R,L,C)

dydt = zeros(2,1);

dydt(1) = ((230*sqrt(2)*cos(100*pi*t))/L)-(R*y(1)/L)-(y(2)/(C*L));
dydt(2) = y(1);


end

--------------------------------------------------------

This is my script:

--------------------------------------------------------

R = ...;
L = ...;
C = ...;

tspan = linspace(0,1,10000);

y0(1) = 0;
y0(2) = 0;

[t,y] = ode45(@(t,y) dv5(t,y,R,L,C),tspan,y0);

plot(t,y(:,1),t,y(:,2))
legend('i','q')

--------------------------------------------------------

Now my question is... Where the hell did I go wrong?

Because the "i" function doesn't like like the derivative of "q" in the slightest.

>Where the hell did I go wrong?
wrote shitty code.
magic numbers
called your function ode5 for some reason
not enough matrices

How do I improve? Please I want to learn.

user don't leave me here. How do I improve my shitty code?

>electrical engineers in charge of doing anything intelligent

Did I really make such a huge mess?

Why are you doing RlC circuits in the time domain? Why not transfer domain?

Math isn't real.
Prove me wrong.

What dude? I'm a fucking newb shitter, please explain in layman terms.

We have to do this project on rlc circuits and our tutor basically gave us some basic equations for a resistor, inductor and capacitor and told us good luck.

What is the circuit you're trying to analyze? Kind of hard to help when all you've given us is a garbled MATLAB code

Resistor, inductor and capacitor, as in pic related.

Also I messed the OP up. I put "sin" everywhere, but it should be "cos".

so what exactly should the solution look like?
also what are your parameters for R,C,R and omega?
i get a similar oscillating solution

I just pick some values for R, C, L. [math] \omega = 100 \pi [/math].

For this particular case I picked R = 20, C = 0.0001, L = 2.

and what exactly is the problem with the numerical solution? what should it look like?

Well so as we know [math] \displaystyle I(t) = \frac{dQ}{dt}(t)[/math], however in the numerical solution I posted in the first post in the thread, it doesn't look like this is the case.

Zoom in on your plot. If the oscillations of q are small enough, your plot could possibly be correct

but it is! zoom in.
other way to varify is to differentiate y1 via finite differences

Yeah the oscillations are very small indeed.

Now that I take a closer look at it, it might even be the derivative and I wasted all of you guys' time.

Anyway, are there better ways to solve ODEs with Matlab?

nah, not really for initial value problems
ode45 is pretty good, but sometimes you might get a stiff problem for which you need an implicit solver like ode15s or ode23s

Thank you, I will have a look at ode15 and ode23 as well.

It's been a while since I've touched matlab but here goes.

don't call your function ode5.
From here it looks like you code won't even work, since the call is for dv5, but your function is ode5

>dydt(1) = ((230*sqrt(2)*cos(100*pi*t))/L)-(R*y(1)/L)-(y(2)/(C*L));
This line doesn't match your equation, but it looks more correct anyway.

Magic numbers are bad coding
Call 230sqrt2 something like A for amplitude

Instead of passing each parameter separately I'd pass it as a matrix (completely a style choice)
param = [L, R, C, A];

nice anonymous function.

tspan should be something like
tspan = [t0, tf];
Let the ode solver choose how to break it up.

These two lines
y0(1) = 0;
y0(2) = 0;
can be made into:
y0 = [0, 0];

>plot(t,y(:,1),t,y(:,2))
Is this even a thing?
Normally I'd do something like
I = y(:,1);
Q = y(:,2);
hold all
plot(t, I)
plot(t, Q)

Check the shape() of your y(:,1) and y(:,2). I always get the wrong direction the first time.

>tspan should be something like
>tspan = [t0, tf];
>Let the ode solver choose how to break it up.
the ode solver will choose how to split up the interval anyways, but it will interpolate the solution if you give him a linspace

There's no reason to interpolate tho...
Unless for some reason the default is really coarse

Plotyy would make that more readable

there is, if the solver has to refine a lot in certain regions, you'll have a few hundred data points to plot basically lying all in the same area.
shit's ineffective man

You can pass ode45 a t-span with intermediate values- the internals of ode45 will still use the set intervals, but the output values will be interpolated with similar to what would have occurred with matlab forcing those subintervals.

Basically OP, read the fucking documentation.

But they're already calculated.
Applying interpolation is the waste.
Unless plot, plots every little point. I figured it'd be optimized not to do that

Try it for yourself- create a struct of ode45's output and see what it uses internally, then try using a set tspan with a vector output and see what it does.

plot plots the whole vector
plotting a million points takes longer than plotting 100 equidistant point.
it's just a small detail, but when you want to see the evolution of a 2/3d pde solution, plotting actually takes some time

Yeah so the "ode5" thing was actually something that got lost in copypasting. I called it that here for convenience.

Tell me why that like of code doesn't match the equation. In the OP I messed it up and accidentally "sin" where I meant "cos".

Also thanks for the tips.

I might be missing something, but
ode5 doesn't match your call dv5?

What do you mean? If you put "ode5" instead of "dv5" it should be good right?