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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1362 - (show annotations)
Mon Dec 17 02:28:16 2007 UTC (11 years, 4 months ago) by gross
File MIME type: text/plain
File size: 4693 byte(s)
and more on FCT solver
1 /* $Id:$ */
2
3 /*******************************************************
4 *
5 * Copyright 2007 by University of Queensland
6 *
7 * http://esscc.uq.edu.au
8 * Primary Business: Queensland, Australia
9 * Licensed under the Open Software License version 3.0
10 * http://www.opensource.org/licenses/osl-3.0.php
11 *
12 *******************************************************/
13
14 /**************************************************************/
15
16 /* Paso: FCTransportProblem */
17
18 /**************************************************************/
19
20 /* Author: l.gross@uq.edu.au */
21
22 /**************************************************************/
23
24
25 #include "Paso.h"
26 #include "SolverFCT.h"
27 #include "PasoUtil.h"
28
29
30 /**************************************************************/
31
32 /* free all memory used by */
33 void Paso_FCTransportProblem_free(Paso_FCTransportProblem* in) {
34 if (in!=NULL) {
35 in->reference_counter--;
36 if (in->reference_counter<=0) {
37
38 Paso_SystemMatrix_free(in->transport_matrix);
39 Paso_SystemMatrix_free(in->flux_matrix);
40 Paso_MPIInfo_free(in->mpi_info);
41
42 MEMFREE(in->lumped_mass_matrix);
43 MEMFREE(in->colorOf);
44 MEMFREE(in->main_iptr);
45 MEMFREE(in);
46 }
47 }
48 }
49
50 Paso_FCTransportProblem* Paso_FCTransportProblem_getReference(Paso_FCTransportProblem* in) {
51 if (in!=NULL) {
52 ++(in->reference_counter);
53 }
54 }
55
56 Paso_SystemMatrix* Paso_FCTransportProblem_borrowTransportMatrix(Paso_FCTransportProblem* in) {
57 return in->transport_matrix;
58 }
59
60 Paso_SystemMatrix* Paso_FCTransportProblem_borrowFluxMatrix(Paso_FCTransportProblem* in) {
61 return in->flux_matrix;
62 }
63
64 double* Paso_FCTransportProblem_borrowLumpedMassMatrix(Paso_FCTransportProblem* in) {
65 return in->lumped_mass_matrix;
66 }
67
68 dim_t Paso_FCTransportProblem_getTotalNumRows(Paso_FCTransportProblem* in) {
69 return Paso_SystemMatrix_getTotalNumRows(in->transport_matrix);
70 }
71
72 Paso_FCTransportProblem* Paso_FCTransportProblem_alloc(double theta, Paso_SystemMatrixPattern *pattern, int block_size
73
74
75 ) {
76 Paso_SystemMatrixType matrix_type=MATRIX_FORMAT_DEFAULT+MATRIX_FORMAT_BLK1; /* at the moment only block size 1 is supported */
77 Paso_FCTransportProblem* out=NULL;
78 dim_t n,i;
79 index_t iptr,iptr_main,k;
80
81 if ((theta<0.) || (theta >1.)) {
82 Paso_setError(TYPE_ERROR,"Paso_FCTransportProblem_alloc: theta needs to be between 0. and. 1.");
83 return NULL;
84 }
85
86 out=MEMALLOC(1,Paso_FCTransportProblem);
87 if (Paso_checkPtr(out)) return NULL;
88
89 out->theta=theta;
90 out->transport_matrix=Paso_SystemMatrix_alloc(matrix_type,pattern,block_size,block_size);
91 Paso_SystemMatrix_allocBuffer(out->transport_matrix);
92 out->flux_matrix=Paso_SystemMatrix_alloc(matrix_type,pattern,block_size,block_size);
93 out->mpi_info=Paso_MPIInfo_getReference(pattern->mpi_info);
94
95 out->colorOf=NULL;
96 out->main_iptr=NULL;
97 out->lumped_mass_matrix=NULL;
98
99 if (Paso_noError()) {
100 n=Paso_SystemMatrix_getTotalNumRows(out->transport_matrix);
101
102 out->colorOf=MEMALLOC(n,index_t);
103 out->main_iptr=MEMALLOC(n,index_t);
104 out->lumped_mass_matrix=MEMALLOC(n,double);
105
106 if ( ! (Paso_checkPtr(out->colorOf) || Paso_checkPtr(out->main_iptr) || Paso_checkPtr(out->lumped_mass_matrix)) ) {
107
108 printf("Paso_SolverFCT_getFCTransportProblem: Revise coloring!!\n");
109 Paso_Pattern_color(pattern->mainPattern,&(out->num_colors),out->colorOf);
110
111
112 /* identify the main diagonals */
113 #pragma omp parallel for schedule(static) private(i,iptr,iptr_main,k)
114 for (i = 0; i < n; ++i) {
115 for (iptr=pattern->mainPattern->ptr[i];iptr<pattern->mainPattern->ptr[i+1]; ++iptr) {
116 iptr_main=pattern->mainPattern->ptr[0]-1;
117 for (iptr=pattern->mainPattern->ptr[i];iptr<pattern->mainPattern->ptr[i+1]; iptr++) {
118 if (pattern->mainPattern->index[iptr]==i) {
119 iptr_main=iptr;
120 break;
121 }
122 }
123 out->main_iptr[i]=iptr_main;
124 if (iptr_main==pattern->mainPattern->ptr[0]-1)
125 Paso_setError(VALUE_ERROR, "Paso_FCTransportProblem_alloc: no main diagonal");
126 }
127 }
128
129 }
130
131 }
132 if (Paso_noError()) {
133 return out;
134 } else {
135 Paso_FCTransportProblem_free(out);
136 return NULL;
137 }
138 }

  ViewVC Help
Powered by ViewVC 1.1.26