/***************************************************************************/
/* Implementation of the mixed boundary condition for a UDS (or multiple): */
/*       q = hC ( phi - PHI_inf )                                          */
/***************************************************************************/

#include "udf.h"
#include "sg.h"  /* needed for the secondary gradient source macro */

/* Names of the user-defined scalar to be used */
enum
  {
    phi1,
    N_REQUIRED_UDS
  };


DEFINE_PROFILE(scalarMixedBC, thread, nv)
{
  /* constants must be specified correctly for the mixed BC */
  real hC, PHI1_inf;

  /* ==================== */

  face_t f;
  
  real A[ND_ND], dG[ND_ND], dr0[ND_ND], es[ND_ND], dr, A_by_es;
  real Af;
  real beta0, gamma;
  real temp1, temp2;
  
  Thread *t0=thread->t0; 
  hC=20.0;
  PHI1_inf=550.0;


  begin_f_loop(f, thread)
    {

     /* identify the cell thread adjacent to the face thread f */
      cell_t c0 = F_C0(f, thread); 
 
      BOUNDARY_FACE_GEOMETRY(f, thread, A, dr, es, A_by_es, dr0);
      Af=NV_MAG(A);
      gamma=C_UDSI_DIFF(c0, t0, phi1);
 
      if (NULLP(T_STORAGE_R_NV(t0, SV_UDSI_G(phi1))))
	beta0=0;  /* if gradient is not allocated and stored yet, bypass
                     the following macro (it happens when case/data files are being read */
      else
        BOUNDARY_SECONDARY_GRADIENT_SOURCE(beta0, SV_UDSI_G(phi1),dG, es, A_by_es, gamma);

  /* temporary variables used in the profile expression */
      temp1=gamma*A_by_es/dr;
      temp2=hC*Af;
      
      
      F_PROFILE(f, thread, nv)
	= (temp1*C_UDSI(c0, t0, phi1)- beta0 + temp2*PHI1_inf)/(temp2 + temp1);
    }
  end_f_loop(f, thread)
  
}
