Hattori Hanzo |
July 14, 2010 06:44 |
Parallelizing UDF
Hi to everyone,
I need to parallelize this UDF that already works in serial mode:
Code:
#include "udf.h"
#define DENS 1056 // kg/m3 (Densità)
#define R 86e8 // Pa/(m3/s) (R(arteriole) + R(venule) + R(vene); by Kim et al 2010)
#define C 3.7e-12 // m3/Pa (Compliance coronarica; by Kim et al 2010)
#define RAP 266.6 // Pa (Right Atrial Pressure; by Ramachandran et al 2009)
#define dt 0.008 // sec (Time Step)
float P_main_old, P_main_new, P_side_old, P_side_new; // pressioni su ramo principale e laterale
float QR_main, QR_side; // portate nelle resistenze
float dP_main, dP_side; // pressioni a cavallo delle capacitanze
float QMMOUT, QMOUT, PMOUT, QMSOUT, QSOUT, PSOUT;
DEFINE_ADJUST(BOUNDARY,domain)
{
FILE *iniziali;
iniziali = fopen("Valori iniziali.txt","r"); // apro il file "Valori iniziali.txt" in modalità "read"
fscanf (iniziali, "%f %f\n%f %f", &P_main_old, &P_main_new, &P_side_old, &P_side_new); // leggo i 4 valori iniziali
fclose (iniziali);
// Calcolo delle portate volumetriche passanti nelle resistenze
QR_main = (P_main_new-RAP)/R;
QR_side = (P_side_new-RAP)/R;
// Calcolo delle pressioni derivanti dalla compliance
dP_main = (QMOUT-QR_main)/C * dt;
dP_side = (QSOUT-QR_side)/C * dt;
// Calcolo della pressione in uscita
P_main_new = P_main_old+dP_main;
P_side_new = P_side_old+dP_side;
P_main_old = P_main_new;
P_side_old = P_side_new;
}
DEFINE_PROFILE(main_out_pressure, thread, nv)
{
face_t f;
QMMOUT=0.0;
QMOUT=0.0;
begin_f_loop (f, thread)
{
QMMOUT += F_FLUX(f, thread); // portata in massa in uscita[kg/s]
}
end_f_loop (f, thread)
QMOUT = QMMOUT/DENS; // portata volumetrica in uscita [m3/s]
PMOUT = P_main_new;
begin_f_loop (f, thread)
{
F_PROFILE(f, thread, nv) = PMOUT;
}
end_f_loop (f, thread)
}
DEFINE_PROFILE(side_out_pressure, thread, nv)
{
face_t f;
QMSOUT=0.0;
QSOUT=0.0;
begin_f_loop (f, thread)
{
QMSOUT += F_FLUX(f, thread); // portata in massa in uscita[kg/s]
}
end_f_loop (f, thread)
QSOUT = QMSOUT/DENS; // portata volumetrica in uscita [m3/s]
PSOUT = P_side_new;
begin_f_loop (f, thread)
{
F_PROFILE(f, thread, nv) = PSOUT;
}
end_f_loop (f, thread)
}
#define INLET_DIAMETER 0.003
DEFINE_PROFILE(inlet_velocity, thread, position)
{
real x[ND_ND];
real coeff,r,v_max,v_mean;
real a0,a1,b1,a2,b2,a3,b3,a4,b4,a5,b5,a6,b6,a7,b7,a8,b8,w; //defining Fourier coefficient
face_t f;
real flow_time = RP_Get_Real ("flow-time");
a0 = 0.2316;
a1 = 0.08323;
b1 = -0.07056;
a2 = 0.08742;
b2 = -0.05437;
a3 = 0.01792;
b3 = -0.004439;
a4 = -0.007452;
b4 = 0.00399;
a5 = 0.001896;
b5 = 0.007417;
a6 = -0.00403;
b6 = -0.004773;
a7 = -0.0000483;
b7 = 0.0008649;
a8 = -0.001566;
b8 = -0.00008407;
w = 7.843;
//defining flow velocity via Fourier series
v_mean=a0 +a1*cos(flow_time*w) + b1*sin(flow_time*w) + a2*cos(2*flow_time*w) + b2*sin(2*flow_time*w) + a3*cos(3*flow_time*w) + b3*sin(3*flow_time*w) + a4*cos(4*flow_time*w) + b4*sin(4*flow_time*w) + a5*cos(5*flow_time*w) + b5*sin(5*flow_time*w) + a6*cos(6*flow_time*w) + b6*sin(6*flow_time*w) + a7*cos(7*flow_time*w) + b7*sin(7*flow_time*w) + a8*cos(8*flow_time*w) + b8*sin(8*flow_time*w);
r = INLET_DIAMETER/2.; //Inlet radius
v_max = 2.*v_mean; //Calculating paraboloid vertex z (max velocity)
coeff = -v_max/pow(r,2.);
begin_f_loop(f, thread)
{
F_CENTROID(x,f,thread);
F_PROFILE(f, thread, position) = coeff*(pow(x[0],2.) + pow(x[1],2)) + v_max;
}
end_f_loop(f, thread)
}
Thanks to everyone can help me.
|