Opened 4 years ago
Last modified 4 years ago
#6090 new defect
model implementation (matlab vs OM)
Reported by: | Owned by: | Mahder Alemseged Gebremedhin | |
---|---|---|---|
Priority: | high | Milestone: | Future |
Component: | Code Generation | Version: | |
Keywords: | Cc: |
Description
Hello,
I implemented an lithium sulfur battery model. First i used matlab "ode15s". I got a solution, but i wanted to extend the model. That's why I now try to implement the model in OpenModelica. If implement the model in OpenModelica I get some errors. It's exactly the running version of matlab, but it's not working in OpenModelica. Does anyone no, why I cannot simulate the model in OpenModelica?
Thank you
Attachments (2)
Change History (14)
by , 4 years ago
Attachment: | Kumaresan_2008.mo added |
---|
comment:1 by , 4 years ago
I cannot test your model since i do not have the EnSysLib
:
Error: Failed to insert class Kumaresan_2008 within EnSysLib.Technologies.Energ yStorages.ElectricalStorage.Lithium_Sulfur_Battery.PhysicalModels;
Could you provide the library or a link to where to get it? Otherwise we can't reproduce the error.
If it is confidential or proprietary you can use our PGP encryption and add it to the ticket so that only OpenModelica developers are able to decrypt it.
You can also just send it to me directly at karim.abdelhak@ fh-bielefeld.de (without space).
comment:2 by , 4 years ago
Also it would be helpful if you provided the exact error messages you are getting from simulating the model.
comment:3 by , 4 years ago
Thank you for the fast reply. You don't need to insert the model to the EnSysLib. Because I don't have reference to my library. You can just delete the "within ...". For simulating the model the input "i_app" should have the value 0.394 for the time intervall. I will attach the error message.
comment:4 by , 4 years ago
Thanks, i can test it now!
The first thing i am noticing is that it is a large nonlinear system that cannot be solved because the determinant of corresponding jacobian is NaN -> the jacobian is singular -> the system is singular. This takes place during the first homotopy step of initialization.
This might have some different origins, for starters i noticed that no index reduction is applied on your model. So the natural question: should there be index reduction from your knowledge?
I will consult some people who designed the nonlinear homotopy initialization and see if we can figure this out!
comment:5 by , 4 years ago
I just tryied to solve the problem and thats why changed some compilation/simulation options...
Thank you again for your help! I'm looking forward to further information.
follow-up: 7 comment:6 by , 4 years ago
maybe the issues is that C
is Real
, because
//17 ...(C[6, 8]/C_ref[6]) ^ (q[6, 6]) = (-80355.1) ^ (-0.5)
comment:7 by , 4 years ago
Replying to anonymous:
maybe the issues is that
C
isReal
, because
//17 ...(C[6, 8]/C_ref[6]) ^ (q[6, 6]) = (-80355.1) ^ (-0.5)
Ah yes that really seems suspicious. Matlab can handle complex numbers, but since we compile everything to C and complex numbers are not a primitive type in C, it cannot handle this case.
Maybe you should use Complex
from the MSL if the values are allowed to be complex.
follow-up: 9 comment:8 by , 4 years ago
Good proposal! But in general this is the concentration of a species inside the battery model. So this value should not be negative. How do you receive the information that (C[6, 8]/C_ref[6]) = -80355.1 ?
As you see in the initial equation block, the start-value is C[i, l] = C_ref[i], so the solution of (C[6, 8]/C_ref[6]) should be 1.
Anyway do i just have to declare C as complex? Or do i have to declare are state variables to complex?
comment:9 by , 4 years ago
Replying to lukas.koenemann@…:
Good proposal! But in general this is the concentration of a species inside the battery model. So this value should not be negative. How do you receive the information that (C[6, 8]/C_ref[6]) = -80355.1 ?
As you see in the initial equation block, the start-value is C[i, l] = C_ref[i], so the solution of (C[6, 8]/C_ref[6]) should be 1.
It seems like they confused the indices, you are right it should always result in 1.
The problem might still occur somewhere else i think. Worth investigating...
Anyway do i just have to declare C as complex? Or do i have to declare are state variables to complex?
In general you have to declare the values that could take complex values at any time (even temporarily during iterative solving processes) as Complex
. You could try and debug your Matlab code to figure out which one needs to be complex. (I would propose to use modified breakpoints that check if the imaginary part of a number is not equal to zero)
If you use the Complex
type the best way is to just try it. If the compiler complains about something, it might be that some other values also have to be complex.
E.g. If you have c*r
where c
is Complex
and r
is Real
you might need to make r
Complex, to have it be consistent in variable types. Although i think this specific case is handled properly.
EDIT:
During simulation the value of C[6, 8]
might change and become negative right? So that could be a future problem anyways.
You could also start by defining everything as Complex
and work your way back from there, that might be easier.
comment:10 by , 4 years ago
Thank you for the reply! I will try to define the variables as complex, but i think i have to change the definition of the parameters too. I will have to check how to define complex matrice. I will inform you about the result.
comment:11 by , 4 years ago
I don't know anything about the details of you model, but a battery model with complex numbers looks very suspicious to me. The only physical models I'm aware of that use complex numbers (besides quantum mechanics stuff) are phasor-based models of AC circuits, but a battery is for sure a DC system.
My understanding is that there may be some problem with your model that causes the argument of the square root to become negative. Matlab does not complain about it and carries on with complex numbers, but that doesn't mean the result makes sense, nor that you should try to mimick this behaviour by using the Complex Modelica package.
I understand C[6,8]
is a concentration, so it should be a Real number between 0 and 1. I don't think making it a complex number is a good idea.
comment:12 by , 4 years ago
I tried to run your model using a recent nightly build of OMC. From what I understand, the initial equations you added give rise to a system fo 134 (!) nonlinear implicit equations, which is then solved by tearing with 10 iteration variables. Some of them are initialized, but some others hava a default zero start value.
proper start-values for some of the following iteration variables might help [1] Real phi_I[4](start=2.39, nominal=1) [10] Real phi_II[3](start=0, nominal=1) [2] Real phi_I[3](start=2.39, nominal=1) [3] Real phi_II[8](start=0, nominal=1) [4] Real phi_II[1](start=0, nominal=1) [5] Real phi_II[7](start=0, nominal=1) [6] Real phi_II[6](start=0, nominal=1) [7] Real phi_II[5](start=0, nominal=1) [8] Real phi_II[4](start=0, nominal=1) [9] Real phi_II[2](start=0, nominal=1)
BTW, it seems that the initialization is failing in a different way than the one you experienced, which is possible if you used an older version of OpenModelica with somewhat different symbolic processing.
I don't know the details of your model in Matlab, but I guess you are explicitly initializing the state variables to their initial values there, so that the initialization problem becomes trivial to solve. I would suggest you to do the same in Modelica, so you avoid these kinds of problems. The initial equations you wrote correspond to a difficult nonlinear problem; you may manage to solve it, but that may require some significant work in setting the right start values or using appropriate simplifications and homotopy (see paper).
We are working on a method to pinpoint bad start values for the Newton solver, see draft paper. The paper is currently under review, the method will eventually be implemented in OpenModelica and could help choosing better start alues, but is unfortunately not yet available.
For concentrations, you may also declare min = 0
and max = 1
, so that additional checks are carried out during the solution process, e.g.
Real C[8, n](each min = 0, each max = 1);
For arrays, you need to use the each
qualifier if you want to apply the scalar value to the entire array.
Kumaresan_2008