SADEL

---

Functions list: manzhuk(5|6|7)


DAE (ODE) solvers:

manzhuk(init, x1, x2, intvls, D)

Arguments:

- init is either a vector of n real initial values, where n is the number of unknowns (or a single scalar initial value, in the case of a single ODE).
- x1 and x2 are real, scalar endpoints of the interval over which the solution to the ODE(s) is evaluated. Initial values in init are the values of the ODE function(s) evaluated at x1.
- intvls is the integer number of discretization intervals used to interpolate the solution function. The number of solution points is the number of intervals + 1.
- D is a vector function of the form D(x,y) specifying the right-hand side of the system

manzhuk(init, x1, x2, intvls, G, J1)
manzhuk(init, x1, x2, intvls, G, J1, J2)


[albumimg]1487[/albumimg] [albumimg]1488[/albumimg] [albumimg]1489[/albumimg] [albumimg]1490[/albumimg]

manzhuk. Kinetic curves.sm (13 КиБ) скачан 88 раз(а).
manzhuk. Van der Pol oscillator 1.sm (10 КиБ) скачан 75 раз(а).
manzhuk. Van der Pol oscillator 2.sm (13 КиБ) скачан 66 раз(а).
manzhuk. Lorenz-Attractor.sm (16 КиБ) скачан 66 раз(а).

manzhuk. Kinetic curves.pdf (38 КиБ) скачан 72 раз(а).
manzhuk. Van der Pol oscillator 1.pdf (46 КиБ) скачан 60 раз(а).
manzhuk. Van der Pol oscillator 2.pdf (62 КиБ) скачан 65 раз(а).
manzhuk. Lorenz-Attractor.pdf (65 КиБ) скачан 66 раз(а).

Manual (rus): DAE-solver-manzhuk.pdf (3 МиБ) скачан 96 раз(а).

Links:

1. SADEL.

Updated. Short form added.

[albumimg]1491[/albumimg] [albumimg]1492[/albumimg] [albumimg]1493[/albumimg] [albumimg]1494[/albumimg] [albumimg]1495[/albumimg]

manzhuk. Kinetic curves (short).sm (11 КиБ) скачан 65 раз(а).
manzhuk. Lorenz-Attractor (short).sm (15 КиБ) скачан 60 раз(а).
manzhuk. Van der Pol oscillator 1 (short).sm (9 КиБ) скачан 58 раз(а).
manzhuk. Van der Pol oscillator 2 (short).sm (12 КиБ) скачан 59 раз(а).
manzhuk. Arenstorf Orbit (short).sm (16 КиБ) скачан 65 раз(а).

manzhuk. Kinetic curves (short).pdf (37 КиБ) скачан 60 раз(а).
manzhuk. Lorenz-Attractor (short).pdf (63 КиБ) скачан 59 раз(а).
manzhuk. Van der Pol oscillator 1 (short).pdf (45 КиБ) скачан 52 раз(а).
manzhuk. Van der Pol oscillator 2 (short).pdf (61 КиБ) скачан 54 раз(а).
manzhuk. Arenstorf Orbit (short).pdf (55 КиБ) скачан 60 раз(а).

Updated. manzhuk() now more stable. Don't use hmx parameter. By default, its value is dt / 2.

Wrote

Updated. Short form added.

Thanks Uni.
The complete Smath code might be educative.
mwode45 solves as well.

Jean

ODE manzhuk Arenstorf Orbit.sm (21 КиБ) скачан 73 раз(а).

Circular/elliptic orbits.

ODE rkfixed Orbit.sm (41 КиБ) скачан 68 раз(а).

Hello uni,

It seems this is a quite interesting ODE-DAE solver. I did not understand what the matrices J1 and J2 are. The manual is quite complicated for me. I just supposed that J1,J2 have something to do with DAE (differential algebraic equations). It seems to me that J1 is the Jacobian matrix regarding the vector of parameters (third argument in the function G()). I suppose J2 is again some kind of Jacobian matrix.

If one of your motivation was to have a DAE solver in SMath, I appreciate that very much.

Looking forward for some more examples from you about using this solver - especially DAE systems.

Regards,
Radovan

I'll add the examples later.

Wrote

If one of your motivation was to have a DAE solver in SMath, I appreciate that very much.

DAE are each specific, I don't think a generalized solver can be designed,
only for odd examples ... I may be wrong , and wish the best.

I have more than 12 examples of how to use it, but, unfortunately, symbolic interpreter slows down the computation of the algorithm. So I decided to show an alternative version of its use, using the Script component. I hope that I will succeed.

Simple test 1.

[albumimg]1496[/albumimg]

manzhuk.test1.sm (11 КиБ) скачан 71 раз(а).
manzhuk.test1.pdf (48 КиБ) скачан 80 раз(а).

Updated.

Changes:
- zero point issue partially fixed;
- added support for the simple D(t,x) definition.

Good afternoon!

I use sadel to calculate the parameters of the transient in the RLC-circuit.
Why doesn't time go by in calculations? file attached model2_sadel_p.sm (25 КиБ) скачан 52 раз(а).

Можно посмотреть на описание задачи в изначальной форме?

Wrote

Можно посмотреть на описание задачи в изначальной форме?

Схема замещения https://i.imgur.com/9p8FrH3.png
Система дифуров, описывающая процессы, происходящие при разряде https://i.imgur.com/Qf2fE8w.png
Взаимная индуктивность https://i.imgur.com/lRdX1Hp.png

Входные параметры электрической цепи рассчитаны, но сомнений не вызывают

Для начала нужно в выражении G(t,x,xp) исправить xp на px. Смущает выражение для второй производной. Кажется, что первая и вторая производные dx/dt должны быть как-то связаны. Я не составлял подобных систем и с ходу не могу сказать что не так. Хорошо бы ещё посмотреть на результат моделирования схемы в специализированной программе, чтобы убедиться, что начальные условия и параметры схемы действительно приводят к какому-то не нулевому результату.

Ещё, мне кажется, что не обязательно использовать библиотеку SADEL для решения этой системы, а можно попробовать использовать любые другие решатели, чтобы избавиться от px-параметра.

Wrote

Для начала нужно в выражении G(t,x,xp) исправить xp на px.

действительно! не заметила. Но конечно это исправление не помогло, программа видимо умнее меня, сама догадалась что имелось ввиду

Wrote

Смущает выражение для второй производной. Кажется, что первая и вторая производные dx/dt должны быть как-то связаны. Я не составлял подобных систем и с ходу не могу сказать что не так.

dx/dt=v, dv/dt=a. Когда составлялась система дифуров вместо второй производной использовалась первая

Wrote

Хорошо бы ещё посмотреть на результат моделирования схемы в специализированной программе, чтобы убедиться, что начальные условия и параметры схемы действительно приводят к какому-то не нулевому результату.

это хорошая идея, попробую. Просто я опираюсь на источники, где описаны именно такая схема, именно такая система дифуров и именно такие входные параметры. Я пытаюсь пока что просто повторить. Это дает определенную уверенность, что должно работать. Хотя не исключено, что где-то упущена ошибка.

То, что в результатах расчета время не движется - это такая реакция программы на то, что она не может посчитать? Не могу понять о чем это говорит

Нужно, чтобы ier был равен нулю, тогда можно делать интерпретации результата. Этого можно достичь, если исправить указанную ошибку и отключить задание минимального шага интегрирования (hmn) или сделать его меньше. Дальше нужно проверить опечатки в системе. Я посмотрю под отладчиком почему время не меняется.

Можно писать вот так, чтобы было нагляднее, но промежуточные вычисления в этом случае замедлят работу функции.

Wrote

Можно писать вот так, чтобы было нагляднее, но промежуточные вычисления в этом случае замедлят работу функции.

вы имеете ввиду нагляднее для того, кто пишет/читает функцию? т.е. другая форма записи?

это здорово - я попробовала так написать, так тоже запускается расчет, правда результат тот же плачевный

По сути это та же форма записи, так как результат работы функции определяется последним выражением блока программы. Определение переменных из вектора k можно вынести из программного блока на лист, т.к. это константы. Ещё нужно всем выражениям указать численный тип оптимизации в свойствах контекстного меню.
Что касается численного моделирования, то можно попробовать посчитать в Proteus данную схемку. Там только нужно подумать как имитировать ключ.

А что за книжка-то?