/[escript]/branches/trilinos_from_5897/dudley/src/Assemble_PDE_Points.cpp
ViewVC logotype

Diff of /branches/trilinos_from_5897/dudley/src/Assemble_PDE_Points.cpp

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 6008 by caltinay, Wed Feb 17 23:53:30 2016 UTC revision 6009 by caltinay, Wed Mar 2 04:13:26 2016 UTC
# Line 15  Line 15 
15  *****************************************************************************/  *****************************************************************************/
16    
17    
18  /************************************************************************************/  /****************************************************************************
19    
20  /*    assembles the system of numEq PDEs into the stiffness matrix S right hand side F  */    Assembles the system of numEqu PDEs into the stiffness matrix S right hand
21  /*    the shape functions for test and solution must be identical */    side F
22    
23          d_dirac_{k,m} u_m and y_dirac_k
24    
25  /*      d_dirac_{k,m} u_m yand _dirac_k */    u has p.numComp components in a 3D domain.
26      The shape functions for test and solution must be identical and
27      row_NS == row_NN.
28    
29  /*    u has p.numComp components in a 3D domain. The shape functions for test and solution must be identical  */    Shape of the coefficients:
 /*    and row_NS == row_NN                                                                                  */  
30    
31  /*    Shape of the coefficients: */        d_dirac = p.numEqu x p.numComp
32          y_dirac = p.numEqu
33    
 /*      d_dirac = p.numEqu x p.numComp  */  
 /*      y_dirac = p.numEqu   */  
34    
35    *****************************************************************************/
 /************************************************************************************/  
   
 #define ESNEEDPYTHON  
 #include "esysUtils/first.h"  
36    
37  #include "Assemble.h"  #include "Assemble.h"
38  #include "Util.h"  #include "Util.h"
 #ifdef _OPENMP  
 #include <omp.h>  
 #endif  
39    
40    namespace dudley {
41    
42  /************************************************************************************/  void Assemble_PDE_Points(const Assemble_Parameters& p,
43                             const Dudley_ElementFile* elements,
44                             escript::ASM_ptr mat, escript::Data& F,
45                             const escript::Data& d_dirac,
46                             const escript::Data& y_dirac)
47    {
48        double* F_p = NULL;
49        if (!F.isEmpty()) {
50            F.requireWrite();
51            F_p = F.getSampleDataRW(0);
52        }
53    
54  void  Dudley_Assemble_PDE_Points(Dudley_Assemble_Parameters p,  #pragma omp parallel
                                  Dudley_ElementFile* elements,  
                                  escript::ASM_ptr mat, escript::Data* F,  
                                  const escript::Data* d_dirac, const escript::Data* y_dirac) {  
   
     index_t color, e, row_index;  
     __const double  *d_dirac_p, *y_dirac_p;  
       
     double *F_p=(requireWrite(F), getSampleDataRW(F,0));    /* use comma, to get around the mixed code and declarations thing */  
   
     #pragma omp parallel private(color, d_dirac_p, y_dirac_p)  
55      {      {
56            for (color=elements->minColor;color<=elements->maxColor;color++) {          for (int color=elements->minColor;color<=elements->maxColor;color++) {
57               /*  open loop over all elements: */              // loop over all elements
58               #pragma omp for private(e) schedule(static)  #pragma omp for
59               for(e=0;e<elements->numElements;e++){              for(index_t e=0; e<elements->numElements; e++) {
60                  if (elements->Color[e]==color) {                  if (elements->Color[e]==color) {
61                        const index_t row_index=p.row_DOF[elements->Nodes[INDEX2(0,e,p.NN)]];
62                        if (!y_dirac.isEmpty()) {
63                            const double* y_dirac_p=y_dirac.getSampleDataRO(e);
64                            Dudley_Util_AddScatter(1, &row_index, p.numEqu,
65                                            y_dirac_p, F_p, p.row_DOF_UpperBound);
66                        }
67                                        
68             d_dirac_p=getSampleDataRO(d_dirac, e);                      if (!d_dirac.isEmpty()) {
69                     y_dirac_p=getSampleDataRO(y_dirac, e);                          const double* d_dirac_p=d_dirac.getSampleDataRO(e);
70                                      Assemble_addToSystemMatrix(mat, 1, &row_index,
71                     row_index=p.row_DOF[elements->Nodes[INDEX2(0,e,p.NN)]];                                 p.numEqu, 1, &row_index, p.numComp, d_dirac_p);
72                                  }
73             if (NULL!=y_dirac_p)  Dudley_Util_AddScatter(1,                  } // end color check
74                                                          &row_index,              } // end element loop
75                                                          p.numEqu,          } // end color loop
76                                                          y_dirac_p,      } // end parallel region
                                                         F_p,  
                                                         p.row_DOF_UpperBound);  
             
                    if (NULL!=d_dirac_p) Dudley_Assemble_addToSystemMatrix(mat,  
                                                                    1,  
                                                                    &row_index,  
                                                                    p.numEqu,  
                                                                    1,  
                                                                    &row_index,  
                                                                    p.numComp,  
                                                                    d_dirac_p);  
                 } /* end color check */  
              } /* end element loop */  
          } /* end color loop */  
    } /* end parallel region */  
77  }  }
78    
79    } // namespace dudley
80    

Legend:
Removed from v.6008  
changed lines
  Added in v.6009

  ViewVC Help
Powered by ViewVC 1.1.26