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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 616 - (show annotations)
Wed Mar 22 02:46:56 2006 UTC (13 years, 7 months ago) by elspeth
File MIME type: text/plain
File size: 6022 byte(s)
Copyright added to more source files.

1 /*
2 ************************************************************
3 * Copyright 2006 by ACcESS MNRF *
4 * *
5 * http://www.access.edu.au *
6 * Primary Business: Queensland, Australia *
7 * Licensed under the Open Software License version 3.0 *
8 * http://www.opensource.org/licenses/osl-3.0.php *
9 * *
10 ************************************************************
11 */
12
13 /**************************************************************/
14
15 /* Updates the element matrices: */
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_{i,j} u_j)_i-(B_{i} u)_i+C_{j} u_j-D u */
22
23 /* u has numComp components. */
24
25 /* Shape of the coefficients: */
26
27 /* A = numDim x numDim */
28 /* B = numDim */
29 /* C = numDim */
30 /* D = 1 */
31
32
33 /**************************************************************/
34
35 /* Author: gross@access.edu.au */
36 /* Vresion: $Id$ */
37
38 /**************************************************************/
39
40 #include "Assemble.h"
41
42 /**************************************************************/
43 void Finley_Assemble_PDEMatrix_Single2(dim_t NS,dim_t numDim,dim_t numQuad,
44 double* S,double* DSDX, double* Vol,
45 dim_t NN, double* EM_S,
46 double* A, bool_t extendedA,
47 double* B, bool_t extendedB,
48 double* C, bool_t extendedC,
49 double* D, bool_t extendedD ) {
50 dim_t s,r,i,j,q;
51 double rtmp;
52
53 for (s=0;s<NS;s++) {
54 for (r=0;r<NS;r++) {
55 /**************************************************************/
56 /* process A: */
57 /**************************************************************/
58 if (NULL!=A) {
59 if (extendedA) {
60 for (i=0;i<numDim;i++) {
61 for (j=0;j<numDim;j++) {
62 for (q=0;q<numQuad;q++) {
63 EM_S[INDEX4(0,0,s,r,1,1,NN)]+=Vol[q]*DSDX[INDEX3(s,i,q,NS,numDim)]*
64 A[INDEX3(i,j,q,numDim,numDim)]*DSDX[INDEX3(r,j,q,NS,numDim)];
65
66 }
67 }
68 }
69 } else {
70 for (i=0;i<numDim;i++) {
71 for (j=0;j<numDim;j++) {
72 if (A[INDEX2(i,j,numDim)]!=0) {
73 rtmp=0;
74 for (q=0;q<numQuad;q++) rtmp+=Vol[q]*DSDX[INDEX3(s,i,q,NS,numDim)]*DSDX[INDEX3(r,j,q,NS,numDim)];
75 EM_S[INDEX4(0,0,s,r,1,1,NN)]+=rtmp*A[INDEX2(i,j,numDim)];
76 }
77 }
78 }
79 }
80 }
81 /**************************************************************/
82 /* process B: */
83 /**************************************************************/
84 if (NULL!=B) {
85 if (extendedB) {
86 for (i=0;i<numDim;i++) {
87 for (q=0;q<numQuad;q++) {
88 EM_S[INDEX4(0,0,s,r,1,1,NN)]+=Vol[q]*DSDX[INDEX3(s,i,q,NS,numDim)]*B[INDEX2(i,q,numDim)]*S[INDEX2(r,q,NS)];
89 }
90 }
91 } else {
92 for (i=0;i<numDim;i++) {
93 if (B[i]!=0) {
94 rtmp=0;
95 for (q=0;q<numQuad;q++) rtmp+=Vol[q]*DSDX[INDEX3(s,i,q,NS,numDim)]*S[INDEX2(r,q,NS)];
96 EM_S[INDEX4(0,0,s,r,1,1,NN)]+=rtmp*B[i];
97 }
98 }
99 }
100 }
101 /**************************************************************/
102 /* process C: */
103 /**************************************************************/
104 if (NULL!=C) {
105 if (extendedC) {
106 for (j=0;j<numDim;j++) {
107 for (q=0;q<numQuad;q++) {
108 EM_S[INDEX4(0,0,s,r,1,1,NN)]+=Vol[q]*S[INDEX2(s,q,NS)]*C[INDEX2(j,q,numDim)]*DSDX[INDEX3(r,j,q,NS,numDim)];
109 }
110 }
111 } else {
112 for (j=0;j<numDim;j++) {
113 if (C[j]!=0) {
114 rtmp=0;
115 for (q=0;q<numQuad;q++) rtmp+=Vol[q]*S[INDEX2(s,q,NS)]*DSDX[INDEX3(r,j,q,NS,numDim)];
116 EM_S[INDEX4(0,0,s,r,1,1,NN)]+=rtmp*C[j];
117 }
118 }
119 }
120 }
121 /************************************************************* */
122 /* process D */
123 /**************************************************************/
124 if (NULL!=D) {
125 if (extendedD) {
126 for (q=0;q<numQuad;q++) {
127 EM_S[INDEX4(0,0,s,r,1,1,NN)]+=Vol[q]*S[INDEX2(s,q,NS)]*D[q]*S[INDEX2(r,q,NS)];
128 }
129 } else {
130 if (D[0]!=0) {
131 rtmp=0;
132 for (q=0;q<numQuad;q++) rtmp+=Vol[q]*S[INDEX2(s,q,NS)]*S[INDEX2(r,q,NS)];
133 EM_S[INDEX4(0,0,s,r,1,1,NN)]+=rtmp*D[0];
134 }
135 }
136 }
137 }
138 }
139 }
140 /*
141 * $Log$
142 * Revision 1.3 2005/09/15 03:44:21 jgs
143 * Merge of development branch dev-02 back to main trunk on 2005-09-15
144 *
145 * Revision 1.2.2.1 2005/09/07 06:26:17 gross
146 * the solver from finley are put into the standalone package paso now
147 *
148 * Revision 1.2 2005/07/08 04:07:46 jgs
149 * Merge of development branch back to main trunk on 2005-07-08
150 *
151 * Revision 1.1.1.1.2.1 2005/06/29 02:34:47 gross
152 * some changes towards 64 integers in finley
153 *
154 * Revision 1.1.1.1 2004/10/26 06:53:57 jgs
155 * initial import of project esys2
156 *
157 * Revision 1.2 2004/07/30 04:37:06 gross
158 * escript and finley are linking now and RecMeshTest.py has been passed
159 *
160 * Revision 1.1 2004/07/02 04:21:13 gross
161 * Finley C code has been included
162 *
163 *
164 */

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26