Plot u vs x for fixed t and y

General Discussion
1

In this runnable code, I want to either plot u versus x for a fixed y and time t, or extract the u vs x data. What should I add to this code to achieve that?
// thermal_convection.edp
//
// Discussion:
//
// Forced + Natural heat convection in a pipe
// Navier−Stokes equations + convection−diffusion on temperature
//
// Licensing:
//
// This code is distributed under the MIT license.
//
// Modified:
//
// 26 June 2015
//
// Author:
//
// Florian De Vuyst
//
// Reference:
//
// Florian De Vuyst,
// Numerical modeling of transport problems using freefem++ software -
// with examples in biology, CFD, traffic flow and energy transfer,
// HAL id: cel-00842234
// Numerical modeling of transport problems using freefem++ software -- with examples in biology, CFD, traffic flow and energy transfer - CEL - Cours en ligne
//
cout << “\n”;
cout << “thermal_convection():\n”;
cout << " FreeFem++ version\n";
cout << " Thermal convection and diffusion by Navier-Stokes flow.\n";
//
// Define the boundary lines.
//
real lx = 0.25;
real Lx = 3.0;
real Ly = 1.0;

border c1 (t=0.0, lx) {x = t; y = 0.0;}
border c2 (t=lx, Lx-lx) {x = t; y = 0.0;}
border c3 (t=Lx-lx,Lx) {x = t; y = 0.0;}
border c4 (t=0.0,Ly) {x = Lx; y = t;}
border c5 (t=Lx,0.0) {x = t; y = Ly;}
border c6 (t=Ly,0.0) {x = 0.0; y = t;}
//
// Define the mesh.
//
mesh Th = buildmesh ( c1(10) + c2(200) + c3(10) + c4(30) + c5(100) + c6(30) );
//
// Plot the mesh.
//
plot ( Th, wait = false, ps = “thermal_convection_mesh.ps” );
//
// Define the finite element spaces:
// UH for velocity fields,
// XH for pressure P, temperature THETA, and specific volume TAU.
//
fespace Uh(Th, P1b );
Uh u, v, uold, vold, uh, vh;
fespace Xh(Th, P1 );
Xh p, ph;
Xh theta, thetaold, thetah;
Xh tau;
//
// Initialize the fields.
//
real uin0 = 0.2;
func uin = uin0 * 4.0 * (y/Ly) * (1.0-y/Ly);
u = uin;
uold = u;
v = 0.0;
vold = v;

real thetain = 20.0;
theta = thetain;
thetaold = theta;
//
// Define the weak form of the heat equation.
//
real dt = 0.05;
real heatflux = 100.0;
real kappa = 0.001;

problem HeatStep ( theta, thetah ) =
int2d (Th) ( theta * thetah / dt )

  • int2d (Th) ( convect ([u,v], -dt, thetaold)*thetah/dt )
  • int2d (Th) ( kappa * dx(theta) * dx(thetah) + kappa * dy(theta) * dy(thetah) )
  • int1d (Th, c2) ( kappa * heatflux * thetah )
  • on ( c6, theta = thetain );
    //
    // Define the weak form of the Navier Stokes equations.
    //
    real gravity = 9.81;
    real nu = 0.001;

problem NavierStokesStep ( [u, v, p], [uh, vh, ph] ) =
int2d (Th) ( u * uh / dt )

  • int2d (Th) ( convect ( [uold, vold], -dt, uold ) * uh/dt )
  • int2d (Th) ( nu * dx(u) * dx(uh) + nu * dy(u) * dy(uh) )
  • int2d (Th) ( tau * dx(p) * uh )
  • int2d (Th) ( v * vh / dt )
  • int2d (Th) ( convect ([uold, vold], -dt, vold) * vh / dt )
  • int2d (Th) ( nu * dx(v) * dx(vh) + nu * dy(v) * dy(vh) )
  • int2d (Th) ( tau * dy(p) * vh )
  • int2d (Th) ( gravity * vh )
  • int2d (Th) ( dx(u) * ph + dy(v) * ph )
  • on ( c6, u=uin, v=0 )
  • on ( c1, c2, c3, c5, u=0, v=0 );
    //
    // Time steps
    //
    real cst = 1.0;
    real dtsnap = 0.5;
    real t = 0.0;
    real ttosnap = 2.0;

for ( int it = 1; it <= 100; it++ )
{
t = it * dt;
cout << it << " T = " << t << “\n”;

HeatStep;

thetaold = theta;
tau = cst * theta;

NavierStokesStep;

uold = u;
vold = v;

plot ( [ u, v ], theta, nbiso = 80, value = false );

if ( ttosnap <= t )
{
cout << " Plotting T = " << t << “\n”;
ttosnap = ttosnap + dtsnap;
plot ( [u,v], theta, nbiso = 80, value = false,
ps = “thermal_convection_”+it+“.ps”,
cmm = "Step “+it );
}
}
//
// Terminate.
//
cout << “\n”;
cout << “thermal_convection():\n”;
cout << " Normal end of execution.\n”;

exit ( 0 );