MODERATION TEAM NOTE: This comment was moved. New location of this comment is at:

http://community.wolfram.com/groups/-/m/t/1470252

Good work. You did a great job, highly appreciated.

Actually, I tired to solve Eq.1 by hand using the same method as mention in this link https://www.math.psu.edu/tseng/class/Math251/Notes-PDE%20pt1.pdf

But I having difficulty to find eigenvalues and eigenfunctions because I have a quadratic term in the equation which gives me a Parabolic Cylinder Functions in the answer. That's why I want to solve it on Mathematica but still confused how to solve?

(1) is a common (ie thermal / slab) PDE that any ODE book discusses. conditions upon it reduce the outer equations solutions.

(2) (3) are a hypothesis of a "systems of PDE". neither ODE book I have shows proofs that prove solutions to systems of PDE (which are more variate) being the same as systems of ODE (which can always be reduce in order). though they solve many slab or wave equations.

While I'm not much studied, I think you would be lucky if any of the PDE solvers gave you all solutions because PDE have many more degrees of freedom in solving than do ODE. You'd need a specific trace on a specific surface for an x,y,z solution: but PDE leave certain trace directions as "variable" and that's bound to trip up solvers ability to show all (relevant) sol'n.

Solvers which are quite complicated code - and Wolfram does show this in the Help you can access: this post suggested by Neil helped me understand the overall "problem" with asking "what can NDSolve do and not do".

Advanced NDSolve tutorial

http://reference.wolfram.com/language/tutorial/NDSolveOverview.html

http://community.wolfram.com/groups/-/m/t/1396110?p_p_auth=lf8ySf9M

I checked: https://www.math.psu.edu/tseng/class/Math251/Notes-PDE%20pt1.pdf, but it has no proof or methods than that of inspection. I understand that http://www.it.uu.se/edu/course/homepage/finmet2/vt14/material/Lab3.pdf modeling of 3space directions upon a surface of gradients (and thus possible anti-derivatives of integrals of these) can be done. It uses computer modeling to show possible directions involved in the parent equation of the partial differential, a backtrace analysis.

There is https://web.math.princeton.edu/~seri/homepage/papers/gws-2006-3.pdf

Wikipedia says Cauchy-Riemann equations and "impact upon imaginary" are used in the form of Laplace's equation - showing even a simple single PDE is a touchy thing. Laplace's equation in cylindrical coordinates looks "odd" and it is: it is formed not by simple translation of coord by conditions of preserving IMZ in such ways as Laplace needed (the handling of singularities needed for his solutions). It requires Cauchy-Riemann considerations and "sacrifices" of some areas to arrive at the form.

https://web.math.princeton.edu/~seri/homepage/papers/gws-2006-3.pdf seems to show principles and perhaps parts of proofs for systems of PDE: but I think in a limited way (a few forms where they can be solved).

http://farside.ph.utexas.edu/teaching/jk1/lectures/Electromagnetism.html shows many uses of PDE

i'd say PDE can be reduced to several ODE, because you can introduce trace equations so each is no longer partial. What i mean is simple. If one PDE is a gradient on a surface any (directional derivative allowed) is possible. It's easy to turn a partial equation into a (set of) ordinary differential equation knowing the original equation. But not knowing the equation(s) being solved for: my head is scratching. Analysis of all possible directions will likely result in finding sets of equations. and I do wonder if some proof in some book shows how that is done.

If anyone can suggest to me a free online source that explains PDE and contains proofs (including proofs for system of PDE reduction), I would appreciate it :).

Can you explain to me the physics of this task? What does this system of equations describe? I'll try to come up with a solution method using Mathematica.

OK! But this is a completely different task than I thought. The first equation is defined in one area, and two others in another area. The solutions of these equations are conjugate for $y = y_p$.

There must be another model. Here we must solve the non-stationary problem, since the heat flux at $y = y_p$ is unknown beforehand. The solution tends to stationary for a long time. The question of the velocity of a fluid - is it the quadratic Poiseuille profile $A1Y^2+A2Y+A3$? The flow rate in a porous medium is constant $U_p$?

It remains for us to find out what to set for the velocity on the boundary $y = 1$ - the condition $u=0$ or is it a free surface?

The problem is not solvable in the combination of boundary conditions that you want to use. I will show an example of a solution of a problem with Neumann conditions on two boundaries y = 0, 1:

A1 = 1; A2 = 1; A3 = 1; A4 = 1; \[Epsilon] = 1; \[Kappa] = 1;
Bi = 1; Subscript[U, p] = 1; eq = {D[Subscript[\[Theta], c][x, y], 
     x] - A4/(\[Epsilon]*(y^2 A1 + y*A2 + A3))*
     D[Subscript[\[Theta], c][x, y], y, y] == 
   NeumannValue[1, y == 1] , 
  D[Subscript[\[Theta], s][x, y], y, y] - 
    Bi*(Subscript[\[Theta], s][x, y] - 
       Subscript[\[Theta], f][x, y]) == 
   0  , \[Kappa]*D[Subscript[\[Theta], f][x, y], y, y] + 
    Bi*(Subscript[\[Theta], s][x, y] - 
       Subscript[\[Theta], f][x, y]) - \[Kappa]*Subscript[U, p]*
     D[Subscript[\[Theta], f][x, y], x] == NeumannValue[1, y == 0]  }; 
ic = {Subscript[\[Theta], s][0, y] == 0, 
   Subscript[\[Theta], f][0, y] == 0, 
   Subscript[\[Theta], c][0, y] == 0};
bc = {DirichletCondition[ {Subscript[\[Theta], s][x, y] == 0, 
     Subscript[\[Theta], f][x, y] == 1} , y == 1], 
   DirichletCondition[ {Subscript[\[Theta], c][x, y] == 0, 
     Subscript[\[Theta], s][x, y] == 0} , y == 0] };   

sol = NDSolve[{eq, ic, bc}, {Subscript[\[Theta], c], 
   Subscript[\[Theta], s], Subscript[\[Theta], f]}, {x, 0, 1}, {y, 0, 
   1}]

{Plot3D[Evaluate[Subscript[\[Theta], c][x, y] /. sol], {x, 0, 1}, {y, 
   0, 1}, Mesh -> None, ColorFunction -> Hue, 
  AxesLabel -> {"x", "y", ""}], 
 Plot3D[Evaluate[Subscript[\[Theta], s][x, y] /. sol], {x, 0, 1}, {y, 
   0, 1}, Mesh -> None, ColorFunction -> Hue, 
  AxesLabel -> {"x", "y", ""}], 
 Plot3D[Evaluate[Subscript[\[Theta], f][x, y] /. sol], {x, 0, 1}, {y, 
   0, 1}, Mesh -> None, ColorFunction -> Hue, 
  AxesLabel -> {"x", "y", ""}]}

OK, I translated into the Wolfram Language and made a code for numerical integration. It is not possible to solve the problem analytically.

A1 = 1; A2 = 1; A3 = 1; A4 = 1; \[Epsilon] = 1; \[Kappa] = 1;
Bi = 1; Subscript[U, p] = 1; eq = {D[Subscript[\[Theta], c][x, y], 
     x] - A4/(\[Epsilon]*(y^2 A1 + y*A2 + A3))*
     D[Subscript[\[Theta], c][x, y], y, y] == 0  , 
  D[Subscript[\[Theta], s][x, y], y, y] - 
    Bi*(Subscript[\[Theta], s][x, y] - 
       Subscript[\[Theta], f][x, y]) == 
   0  , \[Kappa]*D[Subscript[\[Theta], f][x, y], y, y] + 
    Bi*(Subscript[\[Theta], s][x, y] - 
       Subscript[\[Theta], f][x, y]) - \[Kappa]*Subscript[U, p]*
     D[Subscript[\[Theta], f][x, y], x] == 0   }; 
ic = {Subscript[\[Theta], s][0, y] == 0, 
   Subscript[\[Theta], f][0, y] == 0, 
   Subscript[\[Theta], c][0, y] == 0};
bc = {DirichletCondition[ {Subscript[\[Theta], s][x, y] == 1, 
     Subscript[\[Theta], f][x, y] == 0, 
     Subscript[\[Theta], c][x, y] == 0} , y == 1], 
   DirichletCondition[ {Subscript[\[Theta], s][x, y] == 0, 
     Subscript[\[Theta], f][x, y] == 0, 
     Subscript[\[Theta], c][x, y] == 1} , y == 0] };   
sol = NDSolve[{eq, ic, bc}, {Subscript[\[Theta], c], 
   Subscript[\[Theta], s], Subscript[\[Theta], f]}, {x, 0, 1}, {y, 0, 
   1}]

{Plot3D[Evaluate[Subscript[\[Theta], c][x, y] /. sol], {x, 0, 1}, {y, 
   0, 1}, Mesh -> None, ColorFunction -> Hue, 
  AxesLabel -> {"x", "y", ""}], 
 Plot3D[Evaluate[Subscript[\[Theta], s][x, y] /. sol], {x, 0, 1}, {y, 
   0, 1}, Mesh -> None, ColorFunction -> Hue, 
  AxesLabel -> {"x", "y", ""}], 
 Plot3D[Evaluate[Subscript[\[Theta], f][x, y] /. sol], {x, 0, 1}, {y, 
   0, 1}, Mesh -> None, ColorFunction -> Hue, 
  AxesLabel -> {"x", "y", ""}]}

Ok, can you call $Version? I checked the code on 11.01 and 11.3. It works. But on 10.3 and earlier versions it does not work.

Yes, I have. Please see the attached file.

Attachments:

problem.jpg

Do you have the equations in a safer format? I don't open Microsoft formatted files. Sometimes they contain viruses.