/[escript]/trunk/paso/src/SolverFCT.c
ViewVC logotype

Diff of /trunk/paso/src/SolverFCT.c

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

revision 1363 by gross, Mon Dec 17 03:42:03 2007 UTC revision 1370 by gross, Wed Jan 2 09:21:43 2008 UTC
# Line 39  void Paso_FCTransportProblem_free(Paso_F Line 39  void Paso_FCTransportProblem_free(Paso_F
39             Paso_SystemMatrix_free(in->flux_matrix);             Paso_SystemMatrix_free(in->flux_matrix);
40             Paso_MPIInfo_free(in->mpi_info);             Paso_MPIInfo_free(in->mpi_info);
41    
42               MEMFREE(in->u);
43             MEMFREE(in->lumped_mass_matrix);             MEMFREE(in->lumped_mass_matrix);
44             MEMFREE(in->row_sum_flux_matrix);             MEMFREE(in->row_sum_flux_matrix);
45               MEMFREE(in->transport_matrix_diagonal);
46             MEMFREE(in->colorOf);             MEMFREE(in->colorOf);
47             MEMFREE(in->main_iptr);             MEMFREE(in->main_iptr);
48             MEMFREE(in);             MEMFREE(in);
# Line 52  Paso_FCTransportProblem* Paso_FCTranspor Line 54  Paso_FCTransportProblem* Paso_FCTranspor
54       if (in!=NULL) {       if (in!=NULL) {
55          ++(in->reference_counter);          ++(in->reference_counter);
56       }       }
57         return in;
58  }      }    
59    
60  Paso_SystemMatrix* Paso_FCTransportProblem_borrowTransportMatrix(Paso_FCTransportProblem* in) {  Paso_SystemMatrix* Paso_FCTransportProblem_borrowTransportMatrix(Paso_FCTransportProblem* in) {
# Line 70  dim_t Paso_FCTransportProblem_getTotalNu Line 73  dim_t Paso_FCTransportProblem_getTotalNu
73      return Paso_SystemMatrix_getTotalNumRows(in->transport_matrix);      return Paso_SystemMatrix_getTotalNumRows(in->transport_matrix);
74  }  }
75    
76  Paso_FCTransportProblem* Paso_FCTransportProblem_alloc(double theta, Paso_SystemMatrixPattern *pattern, int block_size  Paso_FCTransportProblem* Paso_FCTransportProblem_alloc(double theta, double dt_max, Paso_SystemMatrixPattern *pattern, int block_size
77    
78    
79  ) {  ) {
# Line 88  Paso_FCTransportProblem* Paso_FCTranspor Line 91  Paso_FCTransportProblem* Paso_FCTranspor
91       if (Paso_checkPtr(out)) return NULL;       if (Paso_checkPtr(out)) return NULL;
92    
93       out->theta=theta;       out->theta=theta;
94         out->dt_max=dt_max;
95       out->valid_matrices=FALSE;       out->valid_matrices=FALSE;
96       out->transport_matrix=Paso_SystemMatrix_alloc(matrix_type,pattern,block_size,block_size);       out->transport_matrix=Paso_SystemMatrix_alloc(matrix_type,pattern,block_size,block_size);
97       Paso_SystemMatrix_allocBuffer(out->transport_matrix);       Paso_SystemMatrix_allocBuffer(out->transport_matrix);
# Line 98  Paso_FCTransportProblem* Paso_FCTranspor Line 102  Paso_FCTransportProblem* Paso_FCTranspor
102       out->main_iptr=NULL;       out->main_iptr=NULL;
103       out->lumped_mass_matrix=NULL;       out->lumped_mass_matrix=NULL;
104       out->row_sum_flux_matrix=NULL;       out->row_sum_flux_matrix=NULL;
105         out->transport_matrix_diagonal=NULL;
106    
107       if (Paso_noError()) {       if (Paso_noError()) {
108           n=Paso_SystemMatrix_getTotalNumRows(out->transport_matrix);           n=Paso_SystemMatrix_getTotalNumRows(out->transport_matrix);
# Line 106  Paso_FCTransportProblem* Paso_FCTranspor Line 111  Paso_FCTransportProblem* Paso_FCTranspor
111           out->main_iptr=MEMALLOC(n,index_t);           out->main_iptr=MEMALLOC(n,index_t);
112           out->lumped_mass_matrix=MEMALLOC(n,double);           out->lumped_mass_matrix=MEMALLOC(n,double);
113           out->row_sum_flux_matrix=MEMALLOC(n,double);           out->row_sum_flux_matrix=MEMALLOC(n,double);
114             out->transport_matrix_diagonal=MEMALLOC(n,double);
115             out->u=MEMALLOC(n,double);
116    
117           if ( ! (Paso_checkPtr(out->colorOf) || Paso_checkPtr(out->main_iptr) ||           if ( ! (Paso_checkPtr(out->colorOf) || Paso_checkPtr(out->main_iptr) ||
118                   Paso_checkPtr(out->lumped_mass_matrix) || Paso_checkPtr(out->row_sum_flux_matrix)) ) {                   Paso_checkPtr(out->lumped_mass_matrix) || Paso_checkPtr(out->transport_matrix_diagonal) || Paso_checkPtr(out->row_sum_flux_matrix) || Paso_checkPtr(out->u)) ) {
119                            
120               printf("Paso_SolverFCT_getFCTransportProblem: Revise coloring!!\n");               printf("Paso_SolverFCT_getFCTransportProblem: Revise coloring!!\n");
121               Paso_Pattern_color(pattern->mainPattern,&(out->num_colors),out->colorOf);               Paso_Pattern_color(pattern->mainPattern,&(out->num_colors),out->colorOf);
122    
123                            
124                 #pragma omp parallel for schedule(static) private(i)
125                 for (i = 0; i < n; ++i) {
126                    out->lumped_mass_matrix[i]=0.;
127                    out->row_sum_flux_matrix[i]=0.;
128                    out->u[i]=0.;
129                 }
130               /* identify the main diagonals */               /* identify the main diagonals */
131               #pragma omp parallel for schedule(static) private(i,iptr,iptr_main,k)               #pragma omp parallel for schedule(static) private(i,iptr,iptr_main,k)
132               for (i = 0; i < n; ++i) {               for (i = 0; i < n; ++i) {
                 for (iptr=pattern->mainPattern->ptr[i];iptr<pattern->mainPattern->ptr[i+1]; ++iptr) {  
133                      iptr_main=pattern->mainPattern->ptr[0]-1;                      iptr_main=pattern->mainPattern->ptr[0]-1;
134                      for (iptr=pattern->mainPattern->ptr[i];iptr<pattern->mainPattern->ptr[i+1]; iptr++) {                      for (iptr=pattern->mainPattern->ptr[i];iptr<pattern->mainPattern->ptr[i+1]; iptr++) {
135                            if (pattern->mainPattern->index[iptr]==i) {                            if (pattern->mainPattern->index[iptr]==i) {
# Line 128  Paso_FCTransportProblem* Paso_FCTranspor Line 140  Paso_FCTransportProblem* Paso_FCTranspor
140                      out->main_iptr[i]=iptr_main;                      out->main_iptr[i]=iptr_main;
141                      if (iptr_main==pattern->mainPattern->ptr[0]-1)                      if (iptr_main==pattern->mainPattern->ptr[0]-1)
142                           Paso_setError(VALUE_ERROR, "Paso_FCTransportProblem_alloc: no main diagonal");                           Paso_setError(VALUE_ERROR, "Paso_FCTransportProblem_alloc: no main diagonal");
                 }  
143               }               }
144    
145        }        }
# Line 141  Paso_FCTransportProblem* Paso_FCTranspor Line 152  Paso_FCTransportProblem* Paso_FCTranspor
152       return NULL;       return NULL;
153    }    }
154  }  }
155    
156    void Paso_FCTransportProblem_checkinSolution(Paso_FCTransportProblem* in, double* u) {
157        dim_t i, n;
158      
159        n=Paso_FCTransportProblem_getTotalNumRows(in);
160        #pragma omp parallel for schedule(static) private(i)
161        for (i = 0; i < n; ++i) {
162             in->u[i]=u[i];
163        }
164    }

Legend:
Removed from v.1363  
changed lines
  Added in v.1370

  ViewVC Help
Powered by ViewVC 1.1.26