How to solve system of 3rd order differential equations in Matlab

2 次查看（过去 30 天）
LuisGarcia
LuisGarcia 2018-1-27
评论： LuisGarcia 2018-1-31
I want to solve the following system of differential equations in Matlab for g_a and g_b. I'm using cylindrical coordinates (r, theta) and h, ? and ? are constants. I've read the documentation but I cannot see how I can proceed. If possible, I would like to get an analytical solution - not numerical. Thank you
  6 个评论
显示 4更早的评论
LuisGarcia
LuisGarcia 2018-1-30
@John D'Errico None of the replies were satisfactory and I want to re-formulate the question.
Torsten
Torsten 2018-1-30
编辑：Torsten 2018-1-30
Solve the first ODE for g_a'''(theta) and the second for g_b'''(theta) and then use the usual technique to reduce higher-order ODEs to systems of first-order ODEs:
http://www.math.poly.edu/courses/ma2132/Notes/MA2132EquationsToSystems.pdf
Then use ODE15S to solve.
Best wishes
Torsten.

请先登录，再进行评论。

请先登录，再回答此问题。

采纳的回答

David Goodmanson
David Goodmanson 2018-1-30
Hi Luis,
I am going to assume that r is constant and not a function of theta. One way that could happen is if these equations are modeling small excursions about a fixed r, is that the case?
Anyway, if r is constant, then then an analytic solution is possible. Here is an outline, without doing out the details. Lots of small steps to go through, which is typical. First, in order for the units to work out it looks like both b and mu have units of momentum^2. After multiplying the first equation by (1/h)*r^(-1/2) and the second by (mu/h)*r^(-3/2) and rearranging, you can arrive at
c1*ga' +c2*ga''' = c3*gb +c4*gb''
c5*gb' +c6*gb''' = c7*ga +c8*ga''
where the c's are uncomplicated factors that you can work out. Some of the c's include two dimensionless constants that set the problem,
C = l/(h*r^2) and D = mu/h.
Since all terms include a linear factor of ga or gb, their size is undetermined within an overall multiplicative constant, which could have dimensions. Without including that constant, an analytic solution is
ga = cos(m*theta+phi) gb = A*sin(m*theta+phi) (to solve for m and A)
Phi can be determined later by initial conditions (maybe the overall size can be too). After plugging this in and doing the derivatives of the trig functions, you will get a couple of equations
-c1*m +c2*m^3 = A*(c3 -c4*m^2)
A*(c5*m -c6*m^3) = c7 -c8*m^2
Multiplying the two equations together to eliminate A gives a cubic polynomial in the variable m^2. At this point if not before, you can bring in Matlab to solve this. There will be six possible values for m. For each one you can back substitute into one of the equations to find the appropriate A.
Depending on the constants it could be that m^2 has one positive root which would allow two related physical solutions m = +-sqrt(root) for ga and gb.

更多回答（1 个）

John D'Errico
John D'Errico 2018-1-28
编辑：John D'Errico 2018-1-28
But nobody knows what you have tried. Just saying that you have tried lots of things helps nobody to help you, because that forces someone to write the solution completely for you. In turn, that MIGHT take a great deal of time, because they will need to spend a fair amount of time understanding what you have, what you are trying to solve, what the variables mean, etc.
It appears that you are looking for a solution of this, wanting to find functions g_a(theta) & g_b(theta), for any value of r, as well as the other constants. Clearly there is no solution at r==0, a singularity. And we are given no information about how this process evolves if r were to vary, so effectively we can treat r as a constant, a non-zero one at that since these equations will have no solution at all at the singularity. (You might eventually look at a limiting behavior as r approaches zero though. This is also beyond the scope of Answers, as not being a MATLAB question.)
So one thing you can do is multiply (53) by r^(-1/2). That will simplify it by a fair amount. Or you might try multiplying by just r^(1/2). Again, that will simplify things a bit.
Similarly, do the same for (54). I can't know which is more appropriate, but either might work as well as the other.
I would also be looking for a non-dimensionalization of the problem.
Thus, I have NO idea what h, l, u all mean in context here. What are the units of h, l and u? I might guess that all three have similar units. In fact, I might even guess that these variables have units of length, since they appear in ways where they might also cancel the units of length on r. And since we are told that r is a coordinate variable, it has units of length.
For example, I see the term
(1 + h/u)
Since the constant 1 (one, not l) has no units at all, it makes no physical sense to add h/u to 1, UNLESS h and u have the same units, and therefore h/u is also non-dimensional, thus unit-less itself.
So I would suggest it MIGHT make sense to transform the problem in terms of some fixed ratio for h/u. I might also look at non-dimensionalizing r, so work in terms if r/l or r/h, or even r/u. The correct choice depends up on what all of these variables mean to you, in context.
All of this simplifies the problem. It MIGHT help you to understand some characteristics of the solution, much like understanding how Reynolds number impacts the solution of problems where it applies.
Does any of this apply to how you would solve the problem in MATLAB, or if a solution exists at all? Of course not, beyond possibly making the problem a bit less complex. In fact, nothing about the above comments is at all related to MATLAB. But it will simplify the problem a bit. As well, it will improve your understanding of the problem.
In the end, I might be surprised if an analytical solution exists at all. Or not. That won't happen for me without far more understanding of what you are trying to solve than I am willing to invest the time to gain.
Remember that equations, variables, etc., mean little by themselves. So not until you provide sufficient context and understanding can someone really help you here. And that makes this into something that can only become a long, drawn out consulting issue, but not a question about MATLAB.
Other things you can do are to consider other ways of attacking the problem. For example, can you represent g_a(theta) & g_b(theta) in the form of a series in theta? Think Fourier series, or Chebyshev polynomials for example, or other families of functions that would be appropriate. Again, we are given no context at all, so I am just throwing out some random ideas that MIGHT make sense in the proper context.
The point is, you have much that you can do in the way of mathematical analysis worth investing, BEFORE you just throw this at a differential equation solver and hope a solution magically pops out the end. Again though, that makes this not yet a problem about MATLAB, but a problem of mathematics.

类别

在 Help Center 和 File Exchange 中查找有关 Ordinary Differential Equations 的更多信息

标签

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by