Changeset 7145


Ignore:
Timestamp:
01/23/11 11:43:55 (14 years ago)
Author:
Eric.Larour
Message:

Multi-threaded

Location:
issm/trunk/src/c/modules/InterpFromMesh2dx
Files:
1 added
2 edited

Legend:

Unmodified
Added
Removed

  • issm/trunk/src/c/modules/InterpFromMesh2dx/InterpFromMesh2dx.cpp

    r6964 r7145  
    3030        Vec    vec_incontour=NULL;
    3131        double*    incontour=NULL;
     32
     33        /*threading: */
     34        InterpFromMesh2dxThreadStruct gate;
     35        int num=1;
     36
     37        #ifdef _MULTITHREADING_
     38        num=_NUMTHREADS_;
     39        #endif
    3240
    3341        /*some checks*/
    … …  
    8088        }
    8189
    82         /*Loop over the elements*/
    83         if (debug) printf("      interpolation progress:   %5.2lf %%",0.0);
    84         for (i=0;i<nels_data;i++){
     90        /*initialize thread parameters: */
     91        gate.interpolation_type=interpolation_type;
     92        gate.debug=debug;
     93        gate.nels_data=nels_data;
     94        gate.index_data=index_data;
     95        gate.x_data=x_data;
     96        gate.y_data=y_data;
     97        gate.data=data;
     98        gate.xmin=xmin;
     99        gate.xmax=xmax;
     100        gate.ymin=ymin;
     101        gate.ymax=ymax;
     102        gate.nods_prime=nods_prime;
     103        gate.data_prime=data_prime;
     104        gate.x_prime=x_prime;
     105        gate.y_prime=y_prime;
     106        gate.default_values=default_values;
     107        gate.num_default_values=num_default_values;
     108        gate.incontour=incontour;
    85109
    86                 /*display current iteration*/
    87                 if (debug && fmod((double)i,(double)100)==0) printf("\b\b\b\b\b\b\b%5.2lf %%",(double)i/nels_data*100);
    88 
    89                 /*if there is no point inside the domain, go to next iteration*/
    90                 if ( (x_data[(int)index_data[3*i+0]-1]<xmin) && (x_data[(int)index_data[3*i+1]-1]<xmin) && (x_data[(int)index_data[3*i+2]-1]<xmin)) continue;
    91                 if ( (x_data[(int)index_data[3*i+0]-1]>xmax) && (x_data[(int)index_data[3*i+1]-1]>xmax) && (x_data[(int)index_data[3*i+2]-1]>xmax)) continue;
    92                 if ( (y_data[(int)index_data[3*i+0]-1]<ymin) && (y_data[(int)index_data[3*i+1]-1]<ymin) && (y_data[(int)index_data[3*i+2]-1]<ymin)) continue;
    93                 if ( (y_data[(int)index_data[3*i+0]-1]>ymax) && (y_data[(int)index_data[3*i+1]-1]>ymax) && (y_data[(int)index_data[3*i+2]-1]>ymax)) continue;
    94 
    95                 /*get area of the current element (Jacobian = 2 * area)*/
    96                 //area =x2 * y3 - y2*x3 + x1 * y2 - y1 * x2 + x3 * y1 - y3 * x1;
    97                 area=x_data[(int)index_data[3*i+1]-1]*y_data[(int)index_data[3*i+2]-1]-y_data[(int)index_data[3*i+1]-1]*x_data[(int)index_data[3*i+2]-1]
    98                   +  x_data[(int)index_data[3*i+0]-1]*y_data[(int)index_data[3*i+1]-1]-y_data[(int)index_data[3*i+0]-1]*x_data[(int)index_data[3*i+1]-1]
    99                   +  x_data[(int)index_data[3*i+2]-1]*y_data[(int)index_data[3*i+0]-1]-y_data[(int)index_data[3*i+2]-1]*x_data[(int)index_data[3*i+0]-1];
    100 
    101                 /*loop over the prime nodes*/
    102                 for (j=0;j<nods_prime;j++){
    103 
    104                         if(incontour[j]){
    105 
    106                                 /*Get first area coordinate = det(x-x3  x2-x3 ; y-y3   y2-y3)/area*/
    107                                 area_1=((x_prime[j]-x_data[(int)index_data[3*i+2]-1])*(y_data[(int)index_data[3*i+1]-1]-y_data[(int)index_data[3*i+2]-1])
    108                                                 -  (y_prime[j]-y_data[(int)index_data[3*i+2]-1])*(x_data[(int)index_data[3*i+1]-1]-x_data[(int)index_data[3*i+2]-1]))/area;
    109                                 /*Get second area coordinate =det(x1-x3  x-x3 ; y1-y3   y-y3)/area*/
    110                                 area_2=((x_data[(int)index_data[3*i+0]-1]-x_data[(int)index_data[3*i+2]-1])*(y_prime[j]-y_data[(int)index_data[3*i+2]-1])
    111                                                 - (y_data[(int)index_data[3*i+0]-1]-y_data[(int)index_data[3*i+2]-1])*(x_prime[j]-x_data[(int)index_data[3*i+2]-1]))/area;
    112                                 /*Get third area coordinate = 1-area1-area2*/
    113                                 area_3=1-area_1-area_2;
    114 
    115                                 /*is the current point in the current element?*/
    116                                 if (area_1>=0 && area_2>=0 && area_3>=0){
    117 
    118                                         /*Yes ! compute the value on the point*/
    119                                         if (interpolation_type==1){
    120                                                 /*nodal interpolation*/
    121                                                 data_value=area_1*data[(int)index_data[3*i+0]-1]+area_2*data[(int)index_data[3*i+1]-1]+area_3*data[(int)index_data[3*i+2]-1];
    122                                         }
    123                                         else{
    124                                                 /*element interpolation*/
    125                                                 data_value=data[i];
    126                                         }
    127                                         if (isnan(data_value)){
    128                                                 if(num_default_values==1) data_value=default_values[0];
    129                                                 else data_value=default_values[j];
    130                                         }
    131 
    132                                         /*insert value and go to the next point*/
    133                                         VecSetValue(data_prime,j,data_value,INSERT_VALUES);
    134                                 }
    135                         }
    136                 }
    137         }
    138          if (debug) printf("\b\b\b\b\b\b\b%5.2lf %%\n",100.0);
     110        /*launch the thread manager with InterpFromGridToMeshxt as a core: */
     111        LaunchThread(InterpFromMesh2dxt,(void*)&gate,num);
    139112
    140113        /*Assign output pointers:*/

  • issm/trunk/src/c/modules/InterpFromMesh2dx/InterpFromMesh2dx.h

    r5032 r7145  
    99#include "../../toolkits/toolkits.h"
    1010
     11
     12/*threading: */
     13typedef struct{
     14
     15        int interpolation_type;
     16        bool debug;
     17        int  nels_data;
     18        double* index_data;
     19        double* x_data;
     20        double* y_data;
     21        double* data;
     22        double xmin,xmax;
     23        double ymin,ymax;
     24        int    nods_prime;
     25        Vec    data_prime;
     26        double* x_prime;
     27        double* y_prime;
     28        double* default_values;
     29        int     num_default_values;
     30        double*    incontour;
     31
     32
     33} InterpFromMesh2dxThreadStruct;
     34
    1135int InterpFromMesh2dx( Vec* pdata_prime,double* index_data, double* x_data, double* y_data, int nods_data,int nels_data, double* data, int data_length, double* x_prime, double* y_prime, int nods_prime,
    1236                double* default_values,int num_default_values,Contour** contours,int numcontours);
    1337
     38void* InterpFromMesh2dxt(void* vInterpFromMesh2dxThreadStruct);
     39
    1440#endif /* _INTERPFROMMESH2DX_H */
    1541
Note: See TracChangeset for help on using the changeset viewer.