xref: /phasta/svLS/PRECOND.f (revision 712d3df0b59ebebaaeaea358162c8d2c043c6e08)
1*1e99f302SBen Matthews!     This software is Copyright (c) 2012-2015 The Regents of the
2*1e99f302SBen Matthews!     University of California. All Rights Reserved.
3*1e99f302SBen Matthews!
4*1e99f302SBen Matthews!     Permission to copy and modify this software and its documentation
5*1e99f302SBen Matthews!     for educational, research and non-profit purposes, without fee,
6*1e99f302SBen Matthews!     and without a written agreement is hereby granted, provided that
7*1e99f302SBen Matthews!     the above copyright notice, this paragraph and the following three
8*1e99f302SBen Matthews!     paragraphs appear in all copies.
9*1e99f302SBen Matthews!
10*1e99f302SBen Matthews!     Permission to make commercial use of this software may be obtained
11*1e99f302SBen Matthews!     by contacting:
12*1e99f302SBen Matthews!
13*1e99f302SBen Matthews!     Technology Transfer Office
14*1e99f302SBen Matthews!     9500 Gilman Drive, Mail Code 0910
15*1e99f302SBen Matthews!     University of California
16*1e99f302SBen Matthews!     La Jolla, CA 92093-0910
17*1e99f302SBen Matthews!     (858) 534-5815
18*1e99f302SBen Matthews!     invent@ucsd.edu
19*1e99f302SBen Matthews!
20*1e99f302SBen Matthews!     This software program and documentation are copyrighted by The
21*1e99f302SBen Matthews!     Regents of the University of California. The software program and
22*1e99f302SBen Matthews!     documentation are supplied "as is", without any accompanying
23*1e99f302SBen Matthews!     services from The Regents. The Regents does not warrant that the
24*1e99f302SBen Matthews!     operation of the program will be uninterrupted or error-free. The
25*1e99f302SBen Matthews!     end-user understands that the program was developed for research
26*1e99f302SBen Matthews!     purposes and is advised not to rely exclusively on the program for
27*1e99f302SBen Matthews!     any reason.
28*1e99f302SBen Matthews!
29*1e99f302SBen Matthews!     IN NO EVENT SHALL THE UNIVERSITY OF CALIFORNIA BE LIABLE TO ANY
30*1e99f302SBen Matthews!     PARTY FOR DIRECT, INDIRECT, SPECIAL, INCIDENTAL, OR CONSEQUENTIAL
31*1e99f302SBen Matthews!     DAMAGES, INCLUDING LOST PROFITS, ARISING OUT OF THE USE OF THIS
32*1e99f302SBen Matthews!     SOFTWARE AND ITS DOCUMENTATION, EVEN IF THE UNIVERSITY OF
33*1e99f302SBen Matthews!     CALIFORNIA HAS BEEN ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
34*1e99f302SBen Matthews!     THE UNIVERSITY OF CALIFORNIA SPECIFICALLY DISCLAIMS ANY
35*1e99f302SBen Matthews!     WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES
36*1e99f302SBen Matthews!     OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE. THE
37*1e99f302SBen Matthews!     SOFTWARE PROVIDED HEREUNDER IS ON AN "AS IS" BASIS, AND THE
38*1e99f302SBen Matthews!     UNIVERSITY OF CALIFORNIA HAS NO OBLIGATIONS TO PROVIDE
39*1e99f302SBen Matthews!     MAINTENANCE, SUPPORT, UPDATES, ENHANCEMENTS, OR MODIFICATIONS.
40*1e99f302SBen Matthews
41*1e99f302SBen Matthews      SUBROUTINE PRECOND(nFaces, dof, nNo, nnz, commu, cS, face,
42*1e99f302SBen Matthews     2   rowPtr, colPtr, diagPtr, Val, R, W)
43*1e99f302SBen Matthews
44*1e99f302SBen Matthews      INCLUDE "svLS_STD.h"
45*1e99f302SBen Matthews
46*1e99f302SBen Matthews      INTEGER, INTENT(IN) :: nFaces, dof, nNo, nnz
47*1e99f302SBen Matthews      TYPE(svLS_commuType), INTENT(IN) :: commu
48*1e99f302SBen Matthews      TYPE(svLS_cSType), INTENT(IN) :: cS(commu%nTasks)
49*1e99f302SBen Matthews      TYPE(svLS_faceType), INTENT(INOUT) :: face(nFaces)
50*1e99f302SBen Matthews      INTEGER, INTENT(IN) :: rowPtr(2,nNo), colPtr(nnz), diagPtr(nNo)
51*1e99f302SBen Matthews      REAL*8, INTENT(INOUT) :: Val(dof*dof,nnz), R(dof,nNo)
52*1e99f302SBen Matthews      REAL*8, INTENT(OUT) :: W(dof,nNo)
53*1e99f302SBen Matthews
54*1e99f302SBen Matthews      INTEGER i, j, a, b, d, Ac, faIn
55*1e99f302SBen Matthews
56*1e99f302SBen Matthews      SELECT CASE (dof)
57*1e99f302SBen Matthews      CASE (1)
58*1e99f302SBen Matthews         DO Ac=1, nNo
59*1e99f302SBen Matthews            W(1,Ac) = Val(1,diagPtr(Ac))
60*1e99f302SBen Matthews         END DO
61*1e99f302SBen Matthews      CASE(2)
62*1e99f302SBen Matthews         DO Ac=1, nNo
63*1e99f302SBen Matthews            d       = diagPtr(Ac)
64*1e99f302SBen Matthews            W(1,Ac) = Val(1,d)
65*1e99f302SBen Matthews            W(2,Ac) = Val(4,d)
66*1e99f302SBen Matthews         END DO
67*1e99f302SBen Matthews      CASE(3)
68*1e99f302SBen Matthews         DO Ac=1, nNo
69*1e99f302SBen Matthews            d       = diagPtr(Ac)
70*1e99f302SBen Matthews            W(1,Ac) = Val(1,d)
71*1e99f302SBen Matthews            W(2,Ac) = Val(5,d)
72*1e99f302SBen Matthews            W(3,Ac) = Val(9,d)
73*1e99f302SBen Matthews         END DO
74*1e99f302SBen Matthews      CASE(4)
75*1e99f302SBen Matthews         DO Ac=1, nNo
76*1e99f302SBen Matthews            d       = diagPtr(Ac)
77*1e99f302SBen Matthews            W(1,Ac) = Val(1,d)
78*1e99f302SBen Matthews            W(2,Ac) = Val(6,d)
79*1e99f302SBen Matthews            W(3,Ac) = Val(11,d)
80*1e99f302SBen Matthews            W(4,Ac) = Val(16,d)
81*1e99f302SBen Matthews         END DO
82*1e99f302SBen Matthews      CASE DEFAULT
83*1e99f302SBen Matthews         DO Ac=1, nNo
84*1e99f302SBen Matthews            d = diagPtr(Ac)
85*1e99f302SBen Matthews            DO i=1, dof
86*1e99f302SBen Matthews               W(i,Ac) = Val(i*dof-dof+i,d)
87*1e99f302SBen Matthews            END DO
88*1e99f302SBen Matthews         END DO
89*1e99f302SBen Matthews      END SELECT
90*1e99f302SBen Matthews
91*1e99f302SBen Matthews      CALL COMMUV(dof, nNo, commu, cS, W)
92*1e99f302SBen Matthews
93*1e99f302SBen Matthews      DO Ac=1, nNo
94*1e99f302SBen Matthews         d = diagPtr(Ac)
95*1e99f302SBen Matthews         DO i=1, dof
96*1e99f302SBen Matthews            IF (W(i,Ac) .EQ. 0D0) THEN
97*1e99f302SBen Matthews               W(i,Ac) = 1D0
98*1e99f302SBen Matthews               Val(i*dof-dof+i,d) = 1D0
99*1e99f302SBen Matthews            END IF
100*1e99f302SBen Matthews         END DO
101*1e99f302SBen Matthews      END DO
102*1e99f302SBen Matthews
103*1e99f302SBen Matthews      W = 1D0/SQRT(ABS(W))
104*1e99f302SBen Matthews      DO faIn=1, nFaces
105*1e99f302SBen Matthews         IF (.NOT.face(faIn)%incFlag) CYCLE
106*1e99f302SBen Matthews         i = MIN(face(faIn)%dof,dof)
107*1e99f302SBen Matthews         IF (face(faIn)%bGrp .EQ. BC_TYPE_Dir) THEN
108*1e99f302SBen Matthews            DO a=1, face(faIn)%nNo
109*1e99f302SBen Matthews               Ac = face(faIn)%glob(a)
110*1e99f302SBen Matthews               W(1:i,Ac) = W(1:i,Ac)*face(faIn)%val(1:i,a)
111*1e99f302SBen Matthews            END DO
112*1e99f302SBen Matthews         END IF
113*1e99f302SBen Matthews      END DO
114*1e99f302SBen Matthews
115*1e99f302SBen Matthews      SELECT CASE (dof)
116*1e99f302SBen Matthews      CASE (1)
117*1e99f302SBen Matthews         DO Ac=1, nNo
118*1e99f302SBen Matthews            a          = rowPtr(1,Ac)
119*1e99f302SBen Matthews            b          = rowPtr(2,Ac)
120*1e99f302SBen Matthews            Val(1,a:b) = Val(1,a:b)*W(1,Ac)
121*1e99f302SBen Matthews         END DO
122*1e99f302SBen Matthews      CASE(2)
123*1e99f302SBen Matthews         DO Ac=1, nNo
124*1e99f302SBen Matthews            a            = rowPtr(1,Ac)
125*1e99f302SBen Matthews            b            = rowPtr(2,Ac)
126*1e99f302SBen Matthews            Val(1:2,a:b) = Val(1:2,a:b)*W(1,Ac)
127*1e99f302SBen Matthews            Val(3:4,a:b) = Val(3:4,a:b)*W(2,Ac)
128*1e99f302SBen Matthews         END DO
129*1e99f302SBen Matthews      CASE(3)
130*1e99f302SBen Matthews         DO Ac=1, nNo
131*1e99f302SBen Matthews            a            = rowPtr(1,Ac)
132*1e99f302SBen Matthews            b            = rowPtr(2,Ac)
133*1e99f302SBen Matthews            Val(1:3,a:b) = Val(1:3,a:b)*W(1,Ac)
134*1e99f302SBen Matthews            Val(4:6,a:b) = Val(4:6,a:b)*W(2,Ac)
135*1e99f302SBen Matthews            Val(7:9,a:b) = Val(7:9,a:b)*W(3,Ac)
136*1e99f302SBen Matthews         END DO
137*1e99f302SBen Matthews      CASE(4)
138*1e99f302SBen Matthews         DO Ac=1, nNo
139*1e99f302SBen Matthews            a              = rowPtr(1,Ac)
140*1e99f302SBen Matthews            b              = rowPtr(2,Ac)
141*1e99f302SBen Matthews            Val(1:4,a:b)   = Val(1:4,a:b)*W(1,Ac)
142*1e99f302SBen Matthews            Val(5:8,a:b)   = Val(5:8,a:b)*W(2,Ac)
143*1e99f302SBen Matthews            Val(9:12,a:b)  = Val(9:12,a:b)*W(3,Ac)
144*1e99f302SBen Matthews            Val(13:16,a:b) = Val(13:16,a:b)*W(4,Ac)
145*1e99f302SBen Matthews         END DO
146*1e99f302SBen Matthews      CASE DEFAULT
147*1e99f302SBen Matthews         DO Ac=1, nNo
148*1e99f302SBen Matthews            a = rowPtr(1,Ac)
149*1e99f302SBen Matthews            b = rowPtr(2,Ac)
150*1e99f302SBen Matthews            DO i=1, dof
151*1e99f302SBen Matthews               j = i*dof - dof + 1
152*1e99f302SBen Matthews               Val(j:i*dof,a:b) = Val(j:i*dof,a:b)*W(i,Ac)
153*1e99f302SBen Matthews            END DO
154*1e99f302SBen Matthews         END DO
155*1e99f302SBen Matthews      END SELECT
156*1e99f302SBen Matthews      R = W*R
157*1e99f302SBen Matthews      SELECT CASE (dof)
158*1e99f302SBen Matthews      CASE (1)
159*1e99f302SBen Matthews         DO Ac=1, nNo
160*1e99f302SBen Matthews            DO i=rowPtr(1,Ac), rowPtr(2,Ac)
161*1e99f302SBen Matthews               a = colPtr(i)
162*1e99f302SBen Matthews               Val(1,i) = Val(1,i)*W(1,a)
163*1e99f302SBen Matthews            END DO
164*1e99f302SBen Matthews         END DO
165*1e99f302SBen Matthews      CASE (2)
166*1e99f302SBen Matthews         DO Ac=1, nNo
167*1e99f302SBen Matthews            DO i=rowPtr(1,Ac), rowPtr(2,Ac)
168*1e99f302SBen Matthews               a = colPtr(i)
169*1e99f302SBen Matthews               Val(1:3:2,i) = Val(1:3:2,i)*W(1,a)
170*1e99f302SBen Matthews               Val(2:4:2,i) = Val(2:4:2,i)*W(2,a)
171*1e99f302SBen Matthews            END DO
172*1e99f302SBen Matthews         END DO
173*1e99f302SBen Matthews      CASE (3)
174*1e99f302SBen Matthews         DO Ac=1, nNo
175*1e99f302SBen Matthews            DO i=rowPtr(1,Ac), rowPtr(2,Ac)
176*1e99f302SBen Matthews               a = colPtr(i)
177*1e99f302SBen Matthews               Val(1:7:3,i) = Val(1:7:3,i)*W(1,a)
178*1e99f302SBen Matthews               Val(2:8:3,i) = Val(2:8:3,i)*W(2,a)
179*1e99f302SBen Matthews               Val(3:9:3,i) = Val(3:9:3,i)*W(3,a)
180*1e99f302SBen Matthews            END DO
181*1e99f302SBen Matthews         END DO
182*1e99f302SBen Matthews      CASE (4)
183*1e99f302SBen Matthews         DO Ac=1, nNo
184*1e99f302SBen Matthews            DO i=rowPtr(1,Ac), rowPtr(2,Ac)
185*1e99f302SBen Matthews               a = colPtr(i)
186*1e99f302SBen Matthews               Val(1:13:4,i) = Val(1:13:4,i)*W(1,a)
187*1e99f302SBen Matthews               Val(2:14:4,i) = Val(2:14:4,i)*W(2,a)
188*1e99f302SBen Matthews               Val(3:15:4,i) = Val(3:15:4,i)*W(3,a)
189*1e99f302SBen Matthews               Val(4:16:4,i) = Val(4:16:4,i)*W(4,a)
190*1e99f302SBen Matthews            END DO
191*1e99f302SBen Matthews         END DO
192*1e99f302SBen Matthews      CASE DEFAULT
193*1e99f302SBen Matthews         DO Ac=1, nNo
194*1e99f302SBen Matthews            DO i=rowPtr(1,Ac), rowPtr(2,Ac)
195*1e99f302SBen Matthews               a = colPtr(i)
196*1e99f302SBen Matthews               DO b=1, dof
197*1e99f302SBen Matthews                  j = dof*(dof-1) + b
198*1e99f302SBen Matthews                  Val(b:j:dof,i) = Val(b:j:dof,i)*W(b,a)
199*1e99f302SBen Matthews               END DO
200*1e99f302SBen Matthews            END DO
201*1e99f302SBen Matthews         END DO
202*1e99f302SBen Matthews      END SELECT
203*1e99f302SBen Matthews
204*1e99f302SBen Matthews      DO Ac=1, nNo
205*1e99f302SBen Matthews         i = diagPtr(Ac)
206*1e99f302SBen Matthews         DO a=1, dof
207*1e99f302SBen Matthews            b = (a-1)*dof + a
208*1e99f302SBen Matthews            IF (Val(b,i) .EQ. 0D0) Val(b,i) = 1D0
209*1e99f302SBen Matthews         END DO
210*1e99f302SBen Matthews      END DO
211*1e99f302SBen Matthews
212*1e99f302SBen Matthews      DO faIn=1, nFaces
213*1e99f302SBen Matthews         IF (face(faIn)%coupledFlag) THEN
214*1e99f302SBen Matthews            DO a=1, face(faIn)%nNo
215*1e99f302SBen Matthews               Ac = face(faIn)%glob(a)
216*1e99f302SBen Matthews               DO i=1, MIN(face(faIn)%dof,dof)
217*1e99f302SBen Matthews                  face(faIn)%valM(i,a) = face(faIn)%val(i,a)*W(i,Ac)
218*1e99f302SBen Matthews               END DO
219*1e99f302SBen Matthews            END DO
220*1e99f302SBen Matthews         END IF
221*1e99f302SBen Matthews      END DO
222*1e99f302SBen Matthews
223*1e99f302SBen Matthews      RETURN
224*1e99f302SBen Matthews      END SUBROUTINE PRECOND
225