Differential equations

Question

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.

Answer 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.

Answer 2

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.