Differential equation not working

I wanted to implement the equation below:

I wrote the program as follows. But nothing showing up in the plot.

clc; close all; clear all;
t=0.01:0.001:10;
V=575*sin(2*pi*t);
I=100*sin(2*pi*t);
P0=250;
syms V(t)
syms I(t)
ode = diff(V) == (V*I)*diff(I)-((V/0.01)*(((V*I)/P0)-1));
VSol(t) = dsolve(ode);
fplot(VSol(t),[0 10]);

[Torsten]

Why do you specify I and V and then try to solve a differential equation for V ? If you know that V =575*sin(2*pi*t), you don't need to solve a differential equation.

Further, in order to plot a possible solution for V, you will have to specify an initial condition, e.g. V at t=0.

[Reji G]

What should I modify in the code ?

If V(0)=0.1, then ?

I'm not knowing much of these. That's why came to help center

Answer 1

[Torsten]

I = @(t)100*sin(2*pi*t);;
dIdt = @(t) 2*pi*100*cos(2*pi*t);
V0 = 0.1;
fun = @(t,V,I,dIdt) V/I(t)*dIdt(t)-V/0.01*(V*I(t)/250-1.0);
[T,V] = ode45(@(t,V) fun(t,V,I,dIdt),0.01:0.001:10,V0);
plot(T,V)

[Reji G]

Thanks a lot. The shape of waveform is exact. But a modifications required. I tried it, but ending up with error. I want the frequency as 50Hz so, I =100*sin(100*pi*t); What modification need to be done in the code? I'm a beginner.

[Torsten]

Change I and dIdt accordingly.

[Reji G]

I changed I and dIdT to this.

I = @(t)100*sin(100*pi*t);

dIdt = @(t) 10000*pi*cos(100*pi*t);

But There's no plot in the graph. It's empty.

Getting this error

Warning: Failure at t=1.000000e-02. Unable to meet integration tolerances without reducing the step size below the smallest

value allowed (2.775558e-17) at time t.

[Torsten]

As said in another post of yours, you divide by I which is 0 at all integer multiples of 1/200, thus e.g. at t = 1e-2. The right-hand side of your equation gets undefined at these points for t.

You can see the peaks in the solution of your original equation at t = 0.5, 1, 1.5, 2... The solver was generous to integrate over these singular points. Strictly speaking, the results obtained are wrong and the solver had had to stop integration at t=0.5.

[Reji G]

Thank you for the explanation. But I want to plot the graph. Which is available in "Mayr's Model of the Arc Applied to 50- and 60-Hz Interruption in Low-Voltage Devices". They've plotted it somehow

[Walter Roberson]

Interruptions seldom have continuous second derivatives as required by all Runge-Kutta methods. You need to stop the integration at every discontinuity and then restart on the other side of the discontinuity with adjusted boundary conditions.

[Reji G]

Thank you for the explanation. I'm able to understand the meaning of your explanation. But I do not know what changes need to be done in the code. Kindly help with an executable code instead of explanation.

[Torsten]

@Reji G

Which of the 7 figures do you have in mind and how are they obtained ? All of the figure do not ressemble yours.

https://www.semanticscholar.org/paper/Mayr%27s-Model-of-the-Arc-Applied-to-50-and-60-Hz-in-Wright-Wetter/c1dd8d8e3295490bd5f2b54dabf0a08256584ed9

Since I don't have access to the article, I can't tell.

[Walter Roberson]

No can do, but you do not want to hear the explanation of why it cannot be done.

[Reji G]

Fig 2. I've the paper with me. How can I share it to you ?

[Torsten]

The authors seem to have made the same mistake. They just integrated over the singularity.

It's up to you to decide how to restart with I after the singularity. Don't know if it makes sense regarding the physics, but maybe you can just integrate between two singularities and continue I as periodic.

[Reji G]

Okay. What modification need to be done in the code in order to integrate between two singularities and continue I as periodic? Any help interms of code ?

[Torsten]

We can't since you didn't answer how to cope with the singularities. And we also have to proceed ...

[Reji G]

I don't know this much things. In any way if possible to plot the graph, plz help

[Walter Roberson]

Sorry, that cannot be done. You will need to figure it out yourself as you do not want our explanations.

[Torsten]

Okay. What modification need to be done in the code in order to integrate between two singularities and continue I as periodic? Any help interms of code ?

I = @(t)100*sin(2*pi*t);;
dIdt = @(t) 2*pi*100*cos(2*pi*t);
V0 = 0.1;
fun = @(t,V,I,dIdt) V/I(t)*dIdt(t)-V/0.01*(V*I(t)/250-1.0);
[T,V] = ode45(@(t,V) fun(t,V,I,dIdt),0.01:0.001:0.5,V0);
fun_inter = @(t)interp1(T,V,0.01+mod(t-0.01,0.49));
x = -10:0.01:10;
plot(x,fun_inter(x))

[Walter Roberson]

Note that Torsten's code does not integrate between the singularities, and instead runs one cycle and then after that does interpolation on periodic time.

Nearly any result can be obtained when you numerically integrate across a singularity. In some cases you can sensibly define a Cauchy Principal Value https://en.wikipedia.org/wiki/Cauchy_principal_value but it is not clear to me that you can do that in this particular case.

[Reji G]

Thank you Torsten.

For I = @(t)100*sin(2*pi*t); your first code itself is working fine.

My concern is I = @(t)100*sin(100*pi*t); (in place of 1Hz if I give 50 HZ(2->100)), the program is showing 5 errors. I want a plot for I = @(t)100*sin(100*pi*t).

[Walter Roberson]

f = 50;  %Hz
I = @(t)100*sin(2*pi*f*t);
dIdt = @(t) 2*pi*f*cos(2*pi*t);  %guessing here
V0 = 0.1;
fun = @(t,V,I,dIdt) V/I(t)*dIdt(t)-V/0.01*(V*I(t)/250-1.0);
delta = 1/(10*f);
tspan = delta:delta:1/(2*f)-delta;
[T,V] = ode45(@(t,V) fun(t,V,I,dIdt), tspan, V0);
fun_inter = @(t)interp1(T,V,delta+mod(t-delta, 1/f-delta));
x = 0:delta:10;
plot(x,fun_inter(x))

Heavy aliasing on the drawing.

x = 0:delta:1;
plot(x,fun_inter(x))