[Mosheraa]

All,

I've noticed that several other people have posted questions about modeling flow through elastic walls in the past, but I haven't seen any answers.

I'm trying to model blood flow around a stenosis in a compliant artery, using a 2D axisymmetric model with a thin-wall, linearly elastic fluid wall interaction. From what I've read, Fluent definitely has the capability to model this with dynamic mesh. However, I'm a bit stuck on how to correctly implement a DEFINE_GEOMETRY(name,domain,dynamic_thread,positio n) function.

My approach so far is to use a DEFINE_ADJUST function to read pressure data into either a file or user defined memory. Then, I want to use the DEFINE_GEOMETRY to compute the deflection as r_new = r_old*( 1 + D*( P_wall - P_ref ) ) at each point on wall.

I run into problems however, because the DEFINE_GEOMETRY function has some implicit loop over the dynamic face (which I don't know how to access) but the pressure data is accessed/stored through an explicit face loop.

Does anyone have any advice on how to connect the pressure on a face to the displacement of that face? Has anyone had any particular experience using DEFINE_GEOMETRY?

Any thoughts or suggestions would be awesome.

Thanks,

Aaron

[skullCFD]

That's interesting! I used COMSOL or CFD-ACE+ previously for strongly bidirectionally coupled FSI problems. Could you please update the status of your problem here if you get successful with Fluent? Thanks.

[Mosheraa]

Dear SkullCFD,

I have made a little bit of progress with an elastic-wall boundary by using the DEFINE_GRID_MOTION macro instead of DEFINE_GEOM. The code seems to properly adjust the boundary according to the relation

dr/r = Dw*( P - Pref), Dw = r*(1-vp^2)/(E*h)

However, while the boundary seems to be deforming properly, the rest of my mesh is remaining static =(. Hopefully, this just means I need to read a little more documentation about how to set up a dynamic mesh zone.

Another possible problem my code may have, is that it may have convergence problems if it is simultaneously reading pressure and adjusting the boundary with each iteration. This can be easily fixed by using a DEFINE_ADJUST to read pressure data into F_UDMI to prevent any unsteady solver coupling.

As I've said above, my code still has a few bugs, but here's what I have so far:

#include "udf.h"

/*Define Wall Deformation*/
DEFINE_GRID_MOTION( elastic_wall, domain, dynamic_thread, time, dtime)
{
    /*Model Constants*/

        /*Elasticity Parameters*/
        real E = 400000;                                        /*Young's Modulus (Pa)*/
        real vp = .499;                                            /*Poisson's Ratio*/
        real h = .0055;                                            /*Wall Thinkness (m)*/

        /*Reference Values*/
        real Pref = 0;


    /*Variable Declarations*/

        /*Fluent Data Access*/
        Thread *tf = DT_THREAD(dynamic_thread);
        face_t f;
        Node *v;
        int n;

        /*Model Variables*/
        real P;                                                /*Wall Pressure*/
        real r;                                                /*Radius at wall*/
        real NV_VEC(dr);                                  /*Wall Increment*/
        real Dw;                                            /*Wall Distensibility*/


    /*Set deforming flag on adjacent cell zone */

        SET_DEFORMING_THREAD_FLAG(THREAD_T0(tf));            /*??Need to look up exacty what this is doing??*/


    /*Loop Over Faces*/

    begin_f_loop( f, tf )
    {
        /*Compute Pressure on the face*/
        P = -F_P( f, tf );                                    /*??Minus sign to account for face orientation??*/
        Message ("Face Pressure = %f\n", P);

        /*Loop Over Nodes*/
        f_node_loop(f,tf,n)
        {
            /*Assign Node Value*/
            v = F_NODE(f,tf,n);                                

            /*Check if node needs to be updated*/
            if (NODE_POS_NEED_UPDATE (v))
            {
                /* Indicate that node position has been updated*/
                NODE_POS_UPDATED(v);

                /*Compute Wall Distensibility*/
                r = NODE_Y(v);                                /*??Does NODE_Y(v) access radial positiion in axisymmetric??*/
                Dw = r*(1-vp*vp) / (E*h);

                /*Increment Wall*/
                NV_D( dr, =, 0, r*Dw*( P - Pref ),0 );
                NV_V( NODE_COORD(v), +=, dr);

                /*Debugging*/
                Message ("Radius: %f, Deflection: %f\n", r, dr[1] );

            }/*end if statement*/

        }/*end node loop*/

    }/*end face loop*/
    end_f_loop( f, tf )

}/*end DEFINE_GRID_MOTION*/