[brucelqs]

Hi all, I have encountered a problem. After I implanted the anisotropic conductivity UDF in the fluent fuel cell module, and clicked the "Calculate" button, the program did not run, and the following error message popped up when I exited the software "Error at Node 0: No anisotropic UDS diffusivity option" has been selected for difflayer". I tried to add the anisotropic conductivity of UDS-0 (i.e. electron potential) within the source code and through UDF, but both failed. Please help me, I don't know where the problem is.

Here is the code I use inside the source code (pemfc_user.c):

DEFINE_ANISOTROPIC_DIFFUSIVITY(aniso_econd, c, t, i, dmatrix)
{
  anisotropic_elec_cond(c, t, i, dmatrix);
    
  if (pem_initialized && User_Defined_Non_Iso_Econd)
  {
    int zt = get_fc_zone_type(t->id);

    if (zt == adl_zone || zt == cdl_zone || zt == amp_zone || zt == cmp_zone)
    {
      /*--- In this section, users can redefine the "dmatrix" that was computed above for GDL andd MPL only ---*/
        real por, ip, tp, tri, trp;
        por = C_POR(c, t);
        ip = 1. - (3. * por / (2. + por)) * 0.962 * exp(0.367 * (1. - por)) * pow(1. - por, -0.016);
        tp = 1. - (3. * por / (2. + por)) * 0.962 * exp(0.889 * (1. - por)) * pow(1. - por, -0.00715);
        tri = sqrt(ip);
        trp = sqrt(tp);
        real x[3][3];
        static const real diff[3] = { 23000.,23000.,23000. }; //diff is inherent conductivity
        x[1][0] = x[0][1] = 0.0;
        x[2][0] = x[0][2] = 0.0;
        x[2][1] = x[1][2] = 0.0;
        x[2][2] = x[0][0] = tri;
        x[1][1] = trp;

        dmatrix[0][0] = diff[0] * x[0][0] * x[0][0] 
            + diff[1] * x[1][0] * x[1][0]   
            + diff[2] * x[2][0] * x[2][0];  //XX direction,x[1][0]=0,x[2][0]=0
        dmatrix[1][1] = diff[0] * x[0][1] * x[0][1] 
            + diff[1] * x[1][1] * x[1][1]   
            + diff[2] * x[2][1] * x[2][1];  //YY direction,x[0][1]=0,x[2][1]=0
        dmatrix[1][0] = diff[0] * x[0][1] * x[0][0]
            + diff[1] * x[1][1] * x[1][0]
            + diff[2] * x[2][1] * x[2][0];  //dmatrix[1][0]=0
        dmatrix[0][1] = dmatrix[1][0];  //dmatrix[0][1]=0
        dmatrix[2][2] = diff[0] * x[0][2] * x[0][2] 
            + diff[1] * x[1][2] * x[1][2]   
            + diff[2] * x[2][2] * x[2][2];  //ZZ direction,x[0][2]=0,x[1][2]=0
        dmatrix[0][2] = diff[0] * x[0][0] * x[0][2]
            + diff[1] * x[1][0] * x[1][2]
            + diff[2] * x[2][0] * x[2][2];  //dmatrix[0][2]=0
        dmatrix[2][0] = dmatrix[0][2];  //dmatrix[2][0]=0
        dmatrix[1][2] = diff[0] * x[0][1] * x[0][2]
            + diff[1] * x[1][1] * x[1][2]
            + diff[2] * x[2][1] * x[2][2];  //dmatrix[1][2]=0
        dmatrix[2][1] = dmatrix[1][2];  //dmatrix[2][1]=0
    }
  }

  return;
}

I turned on the "Anisotropic E-cond. in GDL/MPL" option in the module, hooked the UDF in the material of GDL and MPL, I tried using the "orthotropic" option, the anisotropic conductivity option was able to calculate correctly , but it fails after loading the UDF.

Below is a screenshot from the software:
1.jpg

2.jpg

3.png

[brucelqs]

Can any kind person help me? I've been stuck with this for a long time and can't think of a solution.

[AlexanderZ]

1. you have a condition

 if (zt == adl_zone || zt == cdl_zone || zt == amp_zone || zt == cmp_zone)

but have only one option, else statement is not defined

2. I don't think its a good idea to define variables inside the loop

 static const real diff[3] = { 23000.,23000.,23000. }; //diff is inherent conductivity

DEFINE_ANISOTROPIC_DIFFUSIVITY macro is a loop by itself

3. I recommend you to start with code from manual. As It should work

[brucelqs]

Can someone kindly help me? I've been stuck with this problem for a few months, and I can pay a fee if you can help me out, thanks a lot!

[AlexanderZ]

show you code

to ask for paid help you may go here
https://www.cfd-online.com/Forums/cfd-freelancers/

[brucelqs]

Here is my code in PEMFC source code:

DEFINE_ANISOTROPIC_DIFFUSIVITY(aniso_econd, c, t, i, dmatrix)
{
anisotropic_elec_cond(c, t, i, dmatrix);

real por, ip, tp;
por = C_POR(c, t);
ip = -0.8155 * por + 0.9549;
tp = -0.9222 * por + 0.7717;

if (pem_initialized && User_Defined_Non_Iso_Econd)
{
int zt = get_fc_zone_type(t->id);

if (zt == adl_zone || zt == cdl_zone || zt == amp_zone || zt == cmp_zone)
{
/*--- In this section, users can redefine the "dmatrix" that was computed above for GDL andd MPL only ---*/
dmatrix[0][0] = 23000. * ip * ip; //XX direction
dmatrix[1][1] = 23000. * tp * tp; //YY direction
dmatrix[1][0] = 0.0; //dmatrix[1][0]=0
dmatrix[0][1] = dmatrix[1][0]; //dmatrix[0][1]=0
dmatrix[2][2] = 23000. * ip * ip; //ZZ direction
dmatrix[0][2] = 0.0; //dmatrix[0][2]=0
dmatrix[2][0] = dmatrix[0][2]; //dmatrix[2][0]=0
dmatrix[1][2] = 0.0; //dmatrix[1][2]=0
dmatrix[2][1] = dmatrix[1][2]; //dmatrix[2][1]=0
}
else
{
dmatrix[0][0] = 0.0;
dmatrix[0][1] = 0.0;
dmatrix[1][0] = 0.0;
dmatrix[1][1] = 0.0;
dmatrix[1][2] = 0.0;
dmatrix[2][1] = 0.0;
dmatrix[2][2] = 0.0;
dmatrix[0][2] = 0.0;
dmatrix[2][0] = 0.0;
}
}

return;
}

When I use this code and compile it, I still have the same problem as before, and Fluent still cannot run.

[AlexanderZ]

did you compile the code?

code seems to be correct.
porosity could be a problem, for test change C_POR(c, t) fo any constant value

[brucelqs]

I replaced C_POR(c,t) with a constant value, but it still doesn't work. I compiled this code from source and it is part of "pemfc_user.c". Then I tried it with UDF and still got the same error. Here is my code to do it with UDF:

#include "udf.h"

DEFINE_ANISOTROPIC_DIFFUSIVITY(aniso_econd_udf, c, t, i, dmatrix)
{
real por = C_POR(c, t);
real ip = 1. - (3. * por / (2. + por)) * 0.962 * exp(0.367 * (1. - por)) * pow(1. - por, -0.016);
real tp = 1. - (3. * por / (2. + por)) * 0.962 * exp(0.889 * (1. - por)) * pow(1. - por, -0.00715);

dmatrix[0][0] = 23000. * ip; //XX direction
dmatrix[1][1] = 23000. * tp; //YY direction
dmatrix[1][0] = 0.0; //dmatrix[1][0]=0
dmatrix[0][1] = dmatrix[1][0]; //dmatrix[0][1]=0
dmatrix[2][2] = 23000. * ip; //ZZ direction
dmatrix[0][2] = 0.0; //dmatrix[0][2]=0
dmatrix[2][0] = dmatrix[0][2]; //dmatrix[2][0]=0
dmatrix[1][2] = 0.0; //dmatrix[1][2]=0
dmatrix[2][1] = dmatrix[1][2]; //dmatrix[2][1]=0
}

I don't know what else I can do.