CFD Online Discussion Forums - UDF for stream function and heatlines

Hello friends

My problem is a 2-D closed cavity with horizontal temperature gradient.
I'm working on a UDF that returns to me the value of stream function and it's profile as well, so far I have partially succeded on this task, because near the walls the profile and values for the stream funtion are not correct (the velocity and stream function should be zero, but I can't get zero because I'm getting the values for the velocities in the center of the cell), besides that, the values and profiles for the rest of the cavity is correct.

Does anybody know how can I get the correct values for the velocities on the walls?

The equation for stream lines is known as psi = psi(i,j-1) + integral(U*dy)
Here's my UDF for the stream funtion:
#include "udf.h"
DEFINE_ON_DEMAND(stream_lines_mem1)
{
Domain *d;
Thread *t;
cell_t c;
real psi=0.;
real vxdy;
real NV_VEC(n);
real NV_VEC(veloc);
real k;
d=Get_Domain(1);

thread_loop_c(t,d)
{
begin_c_loop(c,t)
{
NV_D(veloc,=,C_U(c,t),C_V(c,t));
NV_D(n,=,0.003333,0.0);
vxdy=NV_DOT(n,veloc);
k=C_K_L(c,t);
psi+=vxdy/k;
C_UDMI(c,t,1)=-psi;
}
end_c_loop(c,t)
}
}

The second UDF i'm working on returns the value of the heatlines and it's profile, this UDF is a little more complicated because it considers an inicial value for a control volume (on the bottom for example), and the values for the other volumes is obtained from a integration, in this case, in the vertical direction, the process should be repeated for the next cells until the end of the domain.

The equation for heat lines is known as dH/dy = rho*cp*U*(T-Tref) - k*dT/dx

Does anybody know how can I do this intregration in one direction and consider a different initial value for every control volume on the bottom?

Here's my UDF for heat lines:
#include "udf.h"
DEFINE_ON_DEMAND(heat_lines_mem0)
{
Domain *d;
Thread *t;
cell_t c;
real rho;
real u;
real T;
real k;
real dTdx;
real dTdy;
real H=0.0;
real term1;
real term2;
real Trho;
real aux1;
real aux2;
real dy;
real grad;
real veloc;
d=Get_Domain(1);

thread_loop_c(t,d)
{
begin_c_loop(c,t)
{
if (NULL != THREAD_STORAGE(t,SV_T_G))
{
k=C_K_L(c,t);
T=C_T(c,t);
rho=C_R(c,t);
Trho=T*rho;
veloc=(Trho*C_U(c,t)/k);
grad=(k*C_T_G(c,t)[0]);
dy=0.003333;
aux1=dy*veloc;
aux2=dy*grad;
H+=aux1-aux2;

C_UDMI(c,t,0)=H;
}

else
{
C_UDMI(c,t,0)=0.;
}
}
end_c_loop(c,t)
}
}