Main Page | Directories | File List | File Members

assemble_rhs_B.h File Reference

Go to the source code of this file.

Functions

void assemble_rhs_B (void)


Function Documentation

void assemble_rhs_B void   ) 
 

Definition at line 31 of file assemble_rhs_B.h.

References b, base_geo, base_ref, CONTACT, coord, dim, EPS, i_front, INSULATOR, lg, NA, ND, NE, Ne, nf, ng, Ng, NO, noeud_geo, OHMIC, poids_K, POISSONUPDATEFLAG, Q, SCHOTTKY, and VIC.

Referenced by main().

00032 {
00033  int i,m,l,k1,ni;
00034  double xii[3][3]; // i=1..dim, j=1..lg
00035 
00036  for(k1=0;k1<3;k1++) for(l=0;l<3;l++) xii[k1][l]=0.;
00037  for(i=0;i<=Ng;i++) b[i]=0.;
00038 
00039 // memset(&xii,0,sizeof(xii));
00040 // memset(&b,0,sizeof(b));
00041 
00042 // Assemblage du second membre. See the french book, pag.327, par.9.3.2
00043 // fill the vector b = [f(i,j) f(i+1,j) ...] (non-smart but needed).
00044 // note: current content of vector b is random numbers (not set to zero)
00045  for(m=0;m<Ne;m++){
00046   memset(&xii,0,sizeof(xii));
00047 //  {int i,j; for(i=0;i<3;i++) for(j=0;j<3;j++) xii[i][j]=0.;}
00048   for(l=0;l<lg;l++){
00049 // Evaluer less coordonees cartesiennes du point de Gauss \xi_l
00050    for(k1=0;k1<dim;k1++){ 
00051      int n; 
00052      for(n=0;n<ng;n++) xii[k1][l]+=coord[k1][noeud_geo[n][m]-1]*base_geo[n][l];
00053    }
00054    for(ni=0;ni<nf;ni++){
00055     int n;
00056     double x1;
00057     double dens;
00058     i=noeud_geo[ni][m]-1;
00059     x1=0.;
00060     dens=0.25*(NE[noeud_geo[0][m]-1]+NE[noeud_geo[1][m]-1]
00061               +NE[noeud_geo[2][m]-1]+NE[noeud_geo[3][m]-1]);
00062     for(n=0;n<nf;n++) x1+=-1.e9*Q*(dens-ND[m]+NA[m])*base_ref[n][l];
00063     if(POISSONUPDATEFLAG==NO) x1=0.;
00064 //    for(n=0;n<nf;n++) x1+=sin(2.*PI*coord[2][noeud_geo[n][m]-1]*1.e6)*base_ref[n][l];
00065 //    printf("%g\n",x1);
00066     x1*=base_ref[ni][l];
00067 // Neumann boundary conditions
00068 // We compute \int_{\Gamma_N}{g v}
00069 // see pag.105, formula (4.7) of the french book
00070     if(i_front[i]==INSULATOR) x1+=0.; // i.e. the gradient is zero...
00071     x1*=poids_K[l][m];
00072     b[i]+=x1;
00073 //    printf("value = %g\n",value);
00074     x1=0.;
00075    }
00076   }
00077  }
00078 // non-homogeneous Dirichlet boundary conditions
00079 // *****
00080  for(m=0;m<Ne;m++)
00081   for(ni=0;ni<nf;ni++){
00082     double value;
00083     i=noeud_geo[ni][m]-1;
00084     if(i_front[i]==OHMIC || i_front[i]==CONTACT ||
00085        i_front[i]==SCHOTTKY) value+=VIC[i]/EPS;
00086     b[i]+=value;
00087     value=0.;
00088   }
00089 }


© sourcejam.com 2005-2008