3

I have 4 differential equations. Two 2nd order ODEs

eqns3 = θ1''[y] + Q1 == 0, θ1[0] == θh, -θ1'[0] == Qh
eqns4 = θ2''[y] + Q2 == 0, θ2[1] == 1, -θ2'[1] == Nc (θ2[1] - 1)

and two 4th order ODEs.

eqns1 = θs''''[y] - Bi (k + 1) θs''[y] - Bi k (wf + ws) == 0
eqns2 = θf''''[y] - Bi (k + 1) θf''[y] - Bi k (wf + ws) == 0

The boundary conditions for the 4th order ODEs

θ2[y2] == θf[y2], 
θf[y1] == θ1[y1], 
θ1'[y1] == ke1 θf'[y1] + k ke1 θs'[y1], 
θ2'[y2] == ke2 θf'[y2] + k ke2 θs'[y2]

θs[y1] == θ1[y1], 
θs[y2] == θ2[y2], 
θ1'[y1] == ke1 θf'[y1] + k ke1 θs'[y1], 
θ2'[y2] == ke2 θf'[y2] + k ke2 θs'[y2]

I guess, because the boundary conditions for the fourth order are coupled, Mathematica has been unable to solve for theta f and θs (they are still blue), which does not allow me to plot the graphs.

Also, I have the general solutions of these equations.

θs[y] = Es y^2 + Fs Cosh[y * Sqrt[Bi (k + 1)]] + K1s y + K2s
θf[y] = Ef y^2 + Ff Cosh[y * Sqrt[Bi (k + 1)]] + K1f y + K2f
θ1[y] = A1 y^2 + B1 y + C1
θ1[y] = A2 y^2 + B2 y + C2

Where Es, Fs, K1s, K2s, Ef, Ff, K1f, K2f, A1, B1, C1, A2, B2, C2 are all unknowns. Any tips on how I can solve this problem.

My code

eqns1 = θs''''[y] - Bi (k + 1) θs''[y] - Bi k (wf + ws) == 0
eqns2 = θf''''[y] - Bi (k + 1) θf''[y] - Bi k (wf + ws) == 0
eqns3 = θ1''[y] + Q1 == 0
eqns4 = θ2''[y] + Q2 == 0

DSolve[
  {eqns2, 
   θ2[y2] == θf[y2], 
   θf[y1] == θ1[y1], 
   θ1'[y1] == ke1 θf'[y1] + k ke1 θs'[y1], 
   θ2'[y2] == ke2 θf'[y2] + k ke2 θs'[y2]}, 
  θf[y], y]

DSolve[{eqns3, θ1[0] == θh, -θ1'[0] == Qh}, θ1[y], y]

DSolve[
  {eqns1, 
   θs[y1] == θ1[y1], 
   θs[y2] == θ2[y2], θ1'[y1] == ke1 θf'[y1] + k ke1 θs'[y1], 
   θ2'[y2] == ke2 θf'[y2] + k ke2 θs'[y2]}, 
  θs[y], y]

DSolve[{eqns4, θ2[1] == 1, -θ2'[1] == Nc (θ2[1] - 1)}, θ2[y], y]

Solving for θ1 and θ2 was okay. But for θs and θf it could not solve it. Mathematica gave me an answers for θf and θs but those answers had θf and θs in them.

m_goldberg
  • 107,779
  • 16
  • 103
  • 257
Dammy Ige
  • 31
  • 3
  • Can you supply the code in a workable form? I.e. what you actually tried yourself. – Feyre Aug 21 '16 at 16:02
  • Yes I can. Is there a way I can attach my code file here? – Dammy Ige Aug 21 '16 at 16:22
  • Why not just paste the code in your question? – Feyre Aug 21 '16 at 16:28
  • Just pasted the code into the question – Dammy Ige Aug 21 '16 at 16:35
  • You can use the {} button to format code. – Feyre Aug 21 '16 at 16:50
  • Yeah, thanks. Got it incase of any other questions I have – Dammy Ige Aug 21 '16 at 17:04
  • You might try something like Solve[Eliminate[ts == \[Theta]s[y] /. sol1, \[Theta]s[y]], ts] /. ts -> \[Theta]s[y]. All solutions are very long though, it'd be a lot easier if you could do this numerically. – Feyre Aug 21 '16 at 17:13
  • How do you mean by numerically? Do you mean using NDsolve? Also, what does ts mean in the "Solve[Eliminate[ts == [Theta]s[y] /. sol1, [Theta]s[y]], ts] /. ts -> [Theta]s[y]" ? – Dammy Ige Aug 21 '16 at 17:22
  • Yes, but you'd have to know the values of Bi,k etc. That is an attempt to find theta s without theta s "in them", like you state in the last sentence of your question. – Feyre Aug 21 '16 at 17:26
  • hmmmm, I have the values for {Bi}, {k} et al. So how do I go about it – Dammy Ige Aug 21 '16 at 17:35
  • The solutions I get do not have θf and θs in them: http://i.stack.imgur.com/wxA0e.png -- Are you sure you're using the replacement rules returned by DSolve[] properly? – Michael E2 Aug 21 '16 at 19:36
  • If you have values, plug them in, and just run NDSolve[] the same way, but you need to know ALL variables. – Feyre Aug 21 '16 at 19:50
  • When I used NDSolve for θ1 and θ2 it gave me a result, However the plot do not match what I was expecting. when I used it for θf and θs, it gave me this result – Dammy Ige Aug 22 '16 at 00:30
  • @Feyre When I used NDSolve for θ1 and θ2 it gave me a result, However the plot did not match what I was expecting. When I used it for θf and θs, it said " Boundary conditions not numerical" Also, in the boundary conditions are the derivatives θf' and θs' are the unknowns – Dammy Ige Aug 22 '16 at 00:36
  • @MichaelE2 I have not been using the replacement rule, I have actually been writing them as functions. What is special about the replacement rule and how can I utilize it in this problem. I also looked at your picture and link, please can you explain more – Dammy Ige Aug 22 '16 at 00:41
  • I think I misunderstood. I thought you meant "Mathematica gave me an answers for θf and θs but those answers had θf and θs in them" *respectively*, not that the θf solution had θs in it and vice versa. See my answer for another interpretation. – Michael E2 Aug 22 '16 at 01:41

2 Answers2

1

I think the problem is that some of the equations are coupled, but the OP has written DSolve commands as if they are independent. Try this:

{θfθsSol} = 
  DSolve[DeleteDuplicates@{eqns2, θ2[y2] == θf[y2],
     θf[y1] == θ1[y1], θ1'[y1] == ke1 θf'[y1] + k ke1 θs'[y1],
     θ2'[y2] == ke2 θf'[y2] + k ke2 θs'[y2], 
     eqns1, θs[y1] == θ1[y1], θs[y2] == θ2[y2],
     θ1'[y1] == ke1 θf'[y1] + k ke1 θs'[y1],
     θ2'[y2] == ke2 θf'[y2] + k ke2 θs'[y2]},
  {θf, θs}, y];

One can delete the duplicates by hand, too, if one wants to.

We can check that the solutions are now free of θf and θs with FreeQ:

FreeQ[{θf[y], θs[y]} /. θfθsSol, θf | θs]
(*  True  *)

Note that ReplaceAll (expr /. θfθsSol) is the standard way to substitute the solution into an expression, which one may observe throughout the documentation for DSolve, especially in the plotting examples.

Michael E2
  • 235,386
  • 17
  • 334
  • 747
  • Oh well, maybe all the equations should be put into one big system, although the 3rd & 4th seem to have perfectly self-contained solutions. They could be plugged back into θfθsSol if desired. – Michael E2 Aug 22 '16 at 01:43
  • the 4th order equations were gotten from two second order ODEs [(1/k) [Theta]f''[y] + Bi ([Theta]s - [Theta]f) + wf == 0 [Theta]s''[y] - Bi ([Theta]s - [Theta]f) + ws == 0] that were coupled. After uncoupling, it became a 4th order ODE. I ran the answer you pasted above and did not get any output. and the FreeQ said true. But what does that mean? – Dammy Ige Aug 22 '16 at 12:44
  • @DammyIge Did you read the (linked) documentation page for FreeQ? True means the expression contains no instances of θf or θs (i.e. it's "free" of those symbols). My code gave no (visible) output because of the semicolon. The output is stored in the variable θfθsSol. It a set of two replacement rules that are very long. Execute θfθsSol in a cell by itself to see the output. – Michael E2 Aug 22 '16 at 12:49
  • Now I understand. Executed θfθsSol and still have unknowns C[5] and C[6] in the whole expression just like @JackLaVigne and what I have been getting. Is there by any chance another method mathematically that I can explore? – Dammy Ige Aug 22 '16 at 13:30
  • @DammyIge The reason you have two free constants of integration C[5] and C[6], as Jack points out, is your equations have more degrees of freedom than constraints (boundary conditions). If there are two boundary conditions that should be imposed (that you forgot about, say), then including them would get rid of them. Otherwise, they are free parameters. You can set each to any number and still get a solution. Maybe that makes no sense in your case, but it does in some cases. What to do has to come from the problem itself and the related science or math. In other words, it's up to you. – Michael E2 Aug 22 '16 at 16:50
1

The systems of equations you have in your question are not sufficient to solve without two generated constants. When you look at the equations for θf and θs it looks like there are four pieces of information for each, eight total.

However since they are coupled and two of the equations are identical, there are really only six pieces of information. As a result you get two generated constants.

The work flow path I took was to solve the two 2nd order ODEs first.

DSolveValue[
  {
   θ1''[y] + Q1 == 0,
   θ1'[0] == -Qh,
   θ1[0] == θh
   },
  θ1[y],
  y]

(* 1/2 (-2 Qh y - Q1 y^2 + 2 θh) *)

Define a function for it

θ1fun[y_] := 1/2 (2 θh - 2 Qh y - Q1 y^2)

and

DSolveValue[
  {
   θ2''[y] + Q2 == 0,
   θ2[1] == 1,
   θ2'[1] == Nc (1 - θ2[1])
   },
  θ2[y],
  y]

(* 1/2 (2 - Q2 + 2 Q2 y - Q2 y^2) *)

and a function for it as well

θ2fun[y_] := 1/2 (2 - Q2 + 2 Q2 y - Q2 y^2)

The two coupled equations should be written in one DsolveValue. Use the functions from the solution of the two second order ODE.

{θfy, θsy} = DSolveValue[
  {
   θf''''[y] - Bi (k + 1) θf''[y] - Bi k (wf + ws) == 0,
   θs''''[y] - Bi (k + 1) θs''[y] - Bi k (wf + ws) == 0,
   θf[y2] == θ2fun[y2],
   θf[y1] == θ1fun[y1],
   θ1fun'[y1] == ke1 θf'[y1] + k ke1 θs'[y1],
   θ2fun'[y2] == ke2 θf'[y2] + k ke2 θs'[y2],
   θs[y1] == θ1fun[y1],
   θs[y2] == θ2fun[y2]
   },
  {θf[y], θs[y]},
  y];

This takes a long time. The output is so long that you can't look at it without pain.

You can check however that each has two generated constants.

Cases[θfy, C[i_Integer], Infinity]

The output here has 226 cases of C[5] or C[6].

The second function has fewer occurrences of generated constants.

Cases[θfyθsy[[2]], C[i_Integer], Infinity]

(* {C[5], C[5], C[5], C[5], C[5], C[5], C[6], C[6], C[6], C[6], 
 C[6], C[6]} *)

See if it is possible to re-formulate your problem so the ambiguity is removed. It might make the output simpler.

Jack LaVigne
  • 14,462
  • 2
  • 25
  • 37
  • FWIW, Variables[{θfy, θsy}] should yield all the degrees of freedom, including the constants of integration -- is that what you mean by "ambiguity"? Tho, you may have guessed right, that the OP won't expect or want them, once it's clear that they're there. – Michael E2 Aug 22 '16 at 09:47
  • Thank you guys for your contribution. you are right that the boundary conditions are six instead of eight cos it is duplicated. This problem is from this paper http://www.sciencedirect.com/science/article/pii/S0735193316300537 . The writer used the the general solutions

    { θs[y] = Es y^2 + Fs Cosh[y * Sqrt[Bi (k + 1)]] + K1s y + K2s θf[y] = Ef y^2 + Ff Cosh[y * Sqrt[Bi (k + 1)]] + K1f y + K2f θ1[y] = A1 y^2 + B1 y + C1 θ1[y] = A2 y^2 + B2 y + C2 }

    and somehow arrived at a set of simultaneous equations and used Solve for the unknowns. But it is not clear what he actually did.

    – Dammy Ige Aug 22 '16 at 12:04
  • @DammyIge What if the C1 and C2 are the constants of integration? (Probably not, but I don't have time to look into it right now.) – Michael E2 Aug 22 '16 at 16:51