Helping to solve complex system of DAE by NeuralPDE. #721
Comments
|
@sdesai1287 this would be a good one to try the alternative strategies on. |
|
Ok I can look into that for sure |
|
Using #790 might be simpler. This would be a good problem to dive into. |
|
Hi, I am working on this issue right now.I have tried a the QuasiRandomTraining training strategy with both 50 and 200 points. I plan on trying stochastic training method next, and then going back to the quadature training strategy with different adaptive loss functions. I am also working on the plots right now (I am new to Julia, so I am learning the language as I go. So, the progress is slow. I hope that's fine). I have a question regarding the plot. Which variable corresponds to the phase angle? Thank you! |
@HuynhTran0301, can you answer this? |
|
@sathvikbhagavan In my model, the delta variables represent the rotor phase angle in the power system dynamics. |
|
I am getting this error when I try to plot the graphs of delta1,delta2 and delta3 with respect to the domain (The command used is: plot(domains, delta1) ): Cannot convert Symbolics.VarDomainPairing to series data for plotting. |
|
From what I understand, delta1, delta2, delta3 are all a function of t and domains is defined as: |
|
Do something like: (You have to do a forward pass) ts = 0.0:0.1:25.0 |> collect
x = phi[1](ts', res.u.depvar.delta1)
plot(ts, x')where |
|
I got somewhere close when I used NeuralPDE.QuasiRandomTraining(150) (I attached the image above). QuasiRandomTraining(50) was wildly off, so I will try QuasiRandomTraining(250) to see if I can get something closer to the actual graph. |
|
50 is definitely too small 😅 , that's looking more reasonable. |
|
Right, I am trying to play with different number of points, ranging from small to large, to see how that affects the output. When the number of points is 500, the output is extremely similar to the one I have above: |
|
But I just realised my x-axis needs to be extended, so I am going to run everything again |
|
After rescaling the axis, it seems to match the original OP's result from NeuralPDE. I have trained setting the number of points as: 250, 500, 1000. |
|
Also, in the original code, the domain is given to be domains = [t ∈ Interval(0.0,25.0)] but the graph is plotted from 0 to 2500. The x-axis is given by the domain right? We are plotting rotor phase angle (delta variables) against t? So, shouldn't it be from 0 to 25? |
|
When I plot the figure I just put the data of delta |
|
I have tried the following approaches: QuasiRandomTraining, StochasticTraining with 50, 500 and 1000 points. I also have tried playing the following adaptive loss functions ( paired together with NeuralPDE.QuadratureTraining() ): MiniMaxAdaptiveLoss(5), MiniMaxAdaptiveLoss(10), MiniMaxAdaptiveLoss(20) and NeuralPDE.GradientScaleAdaptiveLoss(5), NeuralPDE.GradientScaleAdaptiveLoss(20). The closest I would get is to QuasiRandomTraining(250) (the graph which I posted a couple of posts above). I am unable to detect the 2nd hump which we can see in the ODE solvers (Tsit5). Is there an alternative route I can try? |
|
Wait, are you using I think other strategies need to be implemented first, @ChrisRackauckas? |
|
Yeah, well that surely would make it bad. This is a good first issue. It should just be made generic and reuse the sampling code of NNODE, I don't see why it wouldn't. |
|
I am just playing with the different strategies in strategy = NeuralPDE.QuadratureTraining() |
|
What do mean by other strategies need to be implemented first? Do you mean in the source code found in https://github.com/SciML/NeuralPDE.jl/blob/master/src/dae_solve.jl#L75 ? I have enjoyed learning and working on this and since Chris mentioned this is a good first issue, I would like to continue working on this as my first issue. |
|
Yes, the strategies need to be refactored such that it can be used for both |
|
I am trying to use NNODE now. Following through the tutorial (https://docs.sciml.ai/NeuralPDE/stable/tutorials/ode/#Solving-an-ODE-with-NNODE), I made the following changes: Everything else is the same. So, in this case, I am getting the following error:
Has it to do with this: Note that NNODE only supports ODEs which are written in the out-of-place form, i.e. du = f(u,p,t), and not f(du,u,p,t). If not declared out-of-place, then the NNODE will exit with an error (https://docs.sciml.ai/NeuralPDE/stable/manual/dae/). But from I see the in |
|
|
|
I am getting the following error when I am implementing the diffeq interface: MethodError: no method matching transform(::Vector{Chain{@NamedTuple{…}, Nothing}}) This is my code: |
|
Instead of chain =[Lux.Chain(Dense(1,n,Lux.tanh),Dense(n,n,Lux.tanh),Dense(n,n,Lux.tanh),Dense(n,n,Lux.tanh),Dense(n,1)) for _ in 1:15]try chain = Lux.Chain(Dense(1,n,Lux.tanh),Dense(n,n,Lux.tanh),Dense(n,n,Lux.tanh),Dense(n,n,Lux.tanh),Dense(n,15))NNODE doesn't support neural network for each variable separately. I think we can add it as a feature. |
|
Hmm, what do you mean add it as a feature? Should I create a pull request and try to add the code and then submit it? Is this because for NeuralPDE.PhysicsInformedNN we have the following positional argument:
So, essentially, for PhysicsInformedNN, the chain has to be a vector but in NNODE and NNADE it can't? Also, after making the changes you suggested above, I am having the error: It seems like the error is in line 69 and |
Yeah, don't worry about it. I will add it soon.
It looks like there is some symbolic stuff somewhere. Looking at your code, this caught my attention: u[10] = -( PeN[1] - (E1^2*Y1[1,1].re + E1*E2*sin(deltaN*(u[1]-u[2]))*Y1[1,2].im + E1*E2*cos(deltaN*(u[1]-u[2]))*Y1[1,2].re + E1*E3*sin(deltaN*(u[1]-u[3]))*Y1[1,3].im + E1*E3*cos(deltaN*(delta1(t)-delta3(t)))*Y1[1,3].re))I think |
|
Hi, I am sorry but I tried spending the last few days making sure that it was symbolic: But I am getting the same error: I am not sure why, I even changed the variables into their specific numbers: E1 = 1.054, E2 = 1.050, E3 = 1.017, omegaN = 120pi, omega_s = 120pi, deltaN = 1, TMN = [0, 1.62342, 0], PeN = [0, 1.62342, 0], PsvN = [0, 1.62342, 0]. |
|
There was a bug in So with using NeuralPDE, Lux, ModelingToolkit, Optimization, OptimizationOptimisers
import ModelingToolkit: Interval
using CSV
using DataFrames
using Plots
using OrdinaryDiffEq
data = CSV.File("/home/sathvikbhagavan/NeuralPDE.jl/test/3gens.csv");
Y1 = CSV.read("/home/sathvikbhagavan/NeuralPDE.jl/test/Y_des.csv", DataFrame, types=Complex{Float64}); #decreasing
# Input of the system.
E1 = 1.054;
E2 = 1.050;
E3 = 1.017;
omegaN = 120*pi;
omega_s = 120*pi
deltaN = 1;
TMN = [0, 1.62342, 0];
PeN = [0, 1.62342, 0];
PsvN = [0, 1.62342, 0];
H = data["H"]
TCH = data["TCH"]
RD = data["RD"]
TSV = data["TSV"]
function eqs(du, u, p, t)
[(u[4]+ 120*pi - 120*pi)/1 - du[1],
(u[5] + 120*pi - 120*pi)/1 - du[2],
(u[6] + 120*pi - 120*pi)/1 - du[3],
(u[7] + 0 - u[10] - 0)/(2*H[1])*120*pi - du[4],
(u[8] + 1.62342 - u[11] - 1.62342)/(2*H[2])*120*pi - du[5],
(u[9] + 0 - u[12] - 0)/(2*H[3])*120*pi - du[6],
-( 0 - ((1.054)^2*Y1[1,1].re + 1.054*1.050*sin(1*(u[1]-u[2]))*Y1[1,2].im + 1.054*1.050*cos(1*(u[1]-u[2]))*Y1[1,2].re + 1.054*1.017*sin(1*(u[1]-u[3]))*Y1[1,3].im + 1.054*1.017*cos(1*(u[1]-u[3]))*Y1[1,3].re)) - du[10],
-( 1.62342 - ((1.050)^2*Y1[2,2].re + 1.054*1.050*sin(1*(u[2]-u[1]))*Y1[2,1].im + 1.054*1.050*cos(1*(u[2]-u[1]))*Y1[2,1].re + 1.050*1.017*sin(1*(u[2]-u[3]))*Y1[2,3].im + 1.050*1.017*cos(1*(u[2]-u[3]))*Y1[2,3].re)) - du[11],
-( 0 - ((1.017)^2*Y1[3,3].re + 1.017*1.054*sin(1*(u[3]-u[1]))*Y1[3,1].im + 1.017*1.054*cos(1*(u[3]-u[1]))*Y1[3,1].re + 1.017*1.050*sin(1*(u[3]-u[2]))*Y1[3,2].im + 1.017*1.050*cos(1*(u[3]-u[2]))*Y1[3,2].re)) - du[12],
(-u[7] - 0 + u[13] + 0) / TCH[1] - du[7],
(-u[8] - 1.62342 + u[14] + 1.62342) / TCH[2] - du[8],
(-u[9] - 0 + u[15] + 0) / TCH[3] - du[9],
(-u[13] - 0 + 0.70945 + 0.33*(-0.0024) - ((u[4] + 120*pi)/120*pi - 1)/RD[1])/TSV[1] - du[13],
(-u[14] - 1.62342 + 1.62342 + 0.334*(-0.0024) - ((u[5] + 120*pi)/120*pi - 1)/RD[2])/TSV[2] - du[14],
(-u[15] - 0 + 0.84843 + 0.33*(-0.0024) - ((u[6] + 120*pi)/120*pi - 1)/RD[3])/TSV[3] - du[15]]
end
bcs = [0.03957/deltaN, 0.3447/deltaN, 0.23038/deltaN, omega_s - omegaN, omega_s - omegaN, omega_s - omegaN, 0.70945 - TMN[1], 1.62342 - TMN[2], 0.848433 - TMN[3], 0.70945 - PeN[1], 1.62342 - PeN[2], 0.848433 - PeN[3], 0.70945 - PsvN[1], 1.62342 - PsvN[2], 0.848433 - PsvN[3]]
dus = zeros(15)
domains = (0.0,25.0)
prob = DAEProblem(eqs, dus, bcs, domains; differential_vars = [[true for i = 1:9]..., false, false, false, [true for i = 1:3]...])
n = 10
chain = Lux.Chain(Dense(1,n,Lux.tanh),Dense(n,n,Lux.tanh),Dense(n,n,Lux.tanh),Dense(n,n,Lux.tanh),Dense(n,15))
alg = NNDAE(chain, OptimizationOptimisers.Adam(0.1))
sol = solve(prob, alg, verbose = true, maxiters = 2000, dt = 1 / 10.0, saveat = 0.01) |
|
Yup, that works. Thanks! I am trying to use the following command to plot the function: However, I am getting the error: I am not sure where the 15-element Vector{Float64} comes from. |
|
I think each element of |
I am trying to solve the DAE system. From the results, although the loss function is very small, the status of the optimization solver is a success, but the results are not the same with the ODE or DAE numerical solvers.
My code is:
The information after solving is:
The following results from NeuralPDE:
And the results from ODE solvers (Tsit5):
The CSV files that I used on the code are attached below.
3gens.csv
Y_des.csv
The text was updated successfully, but these errors were encountered: