/[escript]/trunk/finley/src/Assemble_PDEMatrix_System2.c
ViewVC logotype

Contents of /trunk/finley/src/Assemble_PDEMatrix_System2.c

Parent Directory Parent Directory | Revision Log Revision Log


Revision 471 - (show annotations)
Fri Jan 27 01:33:02 2006 UTC (13 years, 8 months ago) by jgs
File MIME type: text/plain
File size: 8151 byte(s)
reorganise finley src tree to remove inc dir and src/finley directory

1 /*
2 ******************************************************************************
3 * *
4 * COPYRIGHT ACcESS 2003,2004,2005 - All Rights Reserved *
5 * *
6 * This software is the property of ACcESS. No part of this code *
7 * may be copied in any form or by any means without the expressed written *
8 * consent of ACcESS. Copying, use or modification of this software *
9 * by any unauthorised person is illegal unless that person has a software *
10 * license agreement with ACcESS. *
11 * *
12 ******************************************************************************
13 */
14
15 /**************************************************************/
16
17 /* assembles the system of numEq PDEs into the stiffness matrix S: */
18
19 /* -div(A*grad u)-div(B*u)+C*grad u + D*u */
20
21 /* -(A_{k,i,m,j} u_m,j)_i-(B_{k,i,m} u_m)_i+C_{k,m,j} u_m,j-D_{k,m} u_m */
22
23 /* u has numComp components. */
24
25 /* Shape of the coefficients: */
26
27 /* A = numEqu x numDim x numComp x numDim */
28 /* B = numDim x numEqu x numComp */
29 /* C = numEqu x numDim x numComp */
30 /* D = numEqu x numComp */
31
32
33 /**************************************************************/
34
35 /* Author: gross@access.edu.au */
36 /* Version: $Id$ */
37
38 /**************************************************************/
39
40
41 #include "Assemble.h"
42
43 /**************************************************************/
44
45 void Finley_Assemble_PDEMatrix_System2(dim_t NS,dim_t numDim,dim_t numQuad,dim_t numEqu,dim_t numComp,
46 double* S,double* DSDX, double* Vol,
47 int NN, double* EM_S,
48 double* A, bool_t extendedA,
49 double* B, bool_t extendedB,
50 double* C, bool_t extendedC,
51 double* D, bool_t extendedD ) {
52 dim_t s,r,k,m,i,j,q;
53 double rtmp;
54
55
56 /**************************************************************/
57 /* process A: */
58 /**************************************************************/
59 if (NULL!=A) {
60 if (extendedA) {
61 for (s=0;s<NS;s++) {
62 for (r=0;r<NS;r++) {
63 for (k=0;k<numEqu;k++) {
64 for (m=0;m<numComp;m++) {
65 for (i=0;i<numDim;i++) {
66 for (j=0;j<numDim;j++) {
67 for (q=0;q<numQuad;q++) {
68 EM_S[INDEX4(k,m,s,r,numEqu,numComp,NN)]+=Vol[q]*DSDX[INDEX3(s,i,q,NS,numDim)]*
69 A[INDEX5(k,i,m,j,q,numEqu,numDim,numComp,numDim)]*DSDX[INDEX3(r,j,q,NS,numDim)];
70 }
71 }
72 }
73 }
74 }
75 }
76 }
77 } else {
78 for (s=0;s<NS;s++) {
79 for (r=0;r<NS;r++) {
80 for (i=0;i<numDim;i++) {
81 for (j=0;j<numDim;j++) {
82 rtmp=0;
83 for (q=0;q<numQuad;q++) rtmp+=Vol[q]*DSDX[INDEX3(s,i,q,NS,numDim)]*DSDX[INDEX3(r,j,q,NS,numDim)];
84 for (k=0;k<numEqu;k++) {
85 for (m=0;m<numComp;m++) {
86 EM_S[INDEX4(k,m,s,r,numEqu,numComp,NN)]+=rtmp*A[INDEX4(k,i,m,j,numEqu,numDim,numComp)];
87 }
88 }
89 }
90 }
91 }
92 }
93 }
94 }
95 /**************************************************************/
96 /* process B: */
97 /**************************************************************/
98 if (NULL!=B) {
99 if (extendedB) {
100 for (s=0;s<NS;s++) {
101 for (r=0;r<NS;r++) {
102 for (k=0;k<numEqu;k++) {
103 for (m=0;m<numComp;m++) {
104 for (i=0;i<numDim;i++) {
105 for (q=0;q<numQuad;q++) {
106 EM_S[INDEX4(k,m,s,r,numEqu,numComp,NN)]+=Vol[q]*DSDX[INDEX3(s,i,q,NS,numDim)]*
107 B[INDEX4(k,i,m,q,numEqu,numDim,numComp)]*S[INDEX2(r,q,NS)];
108 }
109 }
110 }
111 }
112 }
113 }
114 } else {
115 for (s=0;s<NS;s++) {
116 for (r=0;r<NS;r++) {
117 for (i=0;i<numDim;i++) {
118 rtmp=0;
119 for (q=0;q<numQuad;q++) rtmp+=Vol[q]*DSDX[INDEX3(s,i,q,NS,numDim)]*S[INDEX2(r,q,NS)];
120 for (k=0;k<numEqu;k++) {
121 for (m=0;m<numComp;m++) {
122 EM_S[INDEX4(k,m,s,r,numEqu,numComp,NN)]+=rtmp*B[INDEX3(k,i,m,numEqu,numDim)];
123 }
124 }
125 }
126 }
127 }
128 }
129 }
130 /**************************************************************/
131 /* process C: */
132 /**************************************************************/
133 if (NULL!=C) {
134 if (extendedC) {
135 for (s=0;s<NS;s++) {
136 for (r=0;r<NS;r++) {
137 for (k=0;k<numEqu;k++) {
138 for (m=0;m<numComp;m++) {
139 for (j=0;j<numDim;j++) {
140 for (q=0;q<numQuad;q++) {
141 EM_S[INDEX4(k,m,s,r,numEqu,numComp,NN)]+=Vol[q]*S[INDEX2(s,q,NS)]*
142 C[INDEX4(k,m,j,q,numEqu,numComp,numDim)]*DSDX[INDEX3(r,j,q,NS,numDim)];
143 }
144 }
145 }
146 }
147 }
148 }
149 } else {
150 for (s=0;s<NS;s++) {
151 for (r=0;r<NS;r++) {
152 for (j=0;j<numDim;j++) {
153 rtmp=0;
154 for (q=0;q<numQuad;q++) rtmp+=Vol[q]*S[INDEX2(s,q,NS)]*DSDX[INDEX3(r,j,q,NS,numDim)];
155 for (k=0;k<numEqu;k++) {
156 for (m=0;m<numComp;m++) {
157 EM_S[INDEX4(k,m,s,r,numEqu,numComp,NN)]+=rtmp*C[INDEX3(k,m,j,numEqu,numComp)];
158 }
159 }
160 }
161 }
162 }
163 }
164 }
165 /************************************************************* */
166 /* process D */
167 /**************************************************************/
168 if (NULL!=D) {
169 if (extendedD) {
170 for (s=0;s<NS;s++) {
171 for (r=0;r<NS;r++) {
172 for (k=0;k<numEqu;k++) {
173 for (m=0;m<numComp;m++) {
174 for (q=0;q<numQuad;q++) {
175 EM_S[INDEX4(k,m,s,r,numEqu,numComp,NN)]+=Vol[q]*S[INDEX2(s,q,NS)]*
176 D[INDEX3(k,m,q,numEqu,numComp)]*S[INDEX2(r,q,NS)];
177 }
178 }
179 }
180 }
181 }
182 } else {
183 for (s=0;s<NS;s++) {
184 for (r=0;r<NS;r++) {
185 rtmp=0;
186 for (q=0;q<numQuad;q++) rtmp+=Vol[q]*S[INDEX2(s,q,NS)]*S[INDEX2(r,q,NS)];
187 for (k=0;k<numEqu;k++) {
188 for (m=0;m<numComp;m++) {
189 EM_S[INDEX4(k,m,s,r,numEqu,numComp,NN)]+=rtmp*D[INDEX2(k,m,numEqu)];
190 }
191 }
192 }
193 }
194 }
195 }
196 }
197
198 /*
199 * $Log$
200 * Revision 1.3 2005/09/15 03:44:21 jgs
201 * Merge of development branch dev-02 back to main trunk on 2005-09-15
202 *
203 * Revision 1.2.2.1 2005/09/07 06:26:17 gross
204 * the solver from finley are put into the standalone package paso now
205 *
206 * Revision 1.2 2005/07/08 04:07:46 jgs
207 * Merge of development branch back to main trunk on 2005-07-08
208 *
209 * Revision 1.1.1.1.2.1 2005/06/29 02:34:47 gross
210 * some changes towards 64 integers in finley
211 *
212 * Revision 1.1.1.1 2004/10/26 06:53:57 jgs
213 * initial import of project esys2
214 *
215 * Revision 1.2 2004/07/30 04:37:06 gross
216 * escript and finley are linking now and RecMeshTest.py has been passed
217 *
218 * Revision 1.1 2004/07/02 04:21:13 gross
219 * Finley C code has been included
220 *
221 *
222 */

Properties

Name Value
svn:eol-style native
svn:keywords Author Date Id Revision

  ViewVC Help
Powered by ViewVC 1.1.26