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

Annotation of /trunk/paso/src/Solver_applyBlockDiagonalMatrix.c

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1556 - (hide annotations)
Mon May 12 00:54:58 2008 UTC (12 years, 3 months ago) by gross
File MIME type: text/plain
File size: 2520 byte(s)
Modification to allow mixed mode execution. 
In order to keep the code portable accross platform all MPI calls within
parallel regions have been moved. 


1 ksteube 1312
2 jgs 150 /* $Id$ */
3    
4 ksteube 1312 /*******************************************************
5     *
6     * Copyright 2003-2007 by ACceSS MNRF
7     * Copyright 2007 by University of Queensland
8     *
9     * http://esscc.uq.edu.au
10     * Primary Business: Queensland, Australia
11     * Licensed under the Open Software License version 3.0
12     * http://www.opensource.org/licenses/osl-3.0.php
13     *
14     *******************************************************/
15 dhawcroft 631
16 jgs 150 /**************************************************************/
17    
18     /* Paso: apply block diagonal matrix D: x=D*b */
19    
20     /* should be called within a parallel region */
21     /* barrier synconization should be performed to make sure that the input vector available */
22    
23     /**************************************************************/
24    
25     /* Copyrights by ACcESS Australia 2003, 2004, 2005 */
26     /* Author: gross@access.edu.au */
27    
28     /**************************************************************/
29    
30 gross 700 #include "Paso.h"
31 jgs 150
32     /**************************************************************/
33    
34    
35     void Paso_Solver_applyBlockDiagonalMatrix(dim_t n_block,dim_t n,double* D,index_t* pivot,double* x,double* b) {
36     dim_t i;
37 gross 495 register dim_t i3,i9;
38     register double b0,b1,b2,D00,D10,D20,D01,D11,D21,D02,D12,D22;
39    
40 jgs 150 if (n_block==1) {
41 gross 1556 #pragma omp parallel for private(i) schedule(static)
42 gross 495 for (i=0;i<n;++i) {
43 jgs 150 x[i]=D[i]*b[i];
44     }
45     } else if (n_block==2) {
46 gross 1556 #pragma omp parallel for private(i,b0,b1,D00,D10,D01,D11,i3,i9) schedule(static)
47 gross 495 for (i=0;i<n;++i) {
48     i3=2*i;
49     i9=4*i;
50     b0=b[i3];
51     b1=b[i3+1];
52     D00=D[i9 ];
53     D10=D[i9+1];
54     D01=D[i9+2];
55     D11=D[i9+3];
56     x[i3 ]=D00*b0+D01*b1;
57     x[i3+1]=D10*b0+D11*b1;
58 jgs 150 }
59     } else if (n_block==3) {
60 gross 1556 #pragma omp parallel for private(i,b0,b1,b2,D00,D10,D20,D01,D11,D21,D02,D12,D22,i3,i9) schedule(static)
61 gross 495 for (i=0;i<n;++i) {
62     i3=3*i;
63     i9=9*i;
64     b0=b[i3];
65     b1=b[i3+1];
66     b2=b[i3+2];
67     D00=D[i9 ];
68     D10=D[i9+1];
69     D20=D[i9+2];
70     D01=D[i9+3];
71     D11=D[i9+4];
72     D21=D[i9+5];
73     D02=D[i9+6];
74     D12=D[i9+7];
75     D22=D[i9+8];
76     x[i3 ]=D00*b0+D01*b1+D02*b2;
77     x[i3+1]=D10*b0+D11*b1+D12*b2;
78     x[i3+2]=D20*b0+D21*b1+D22*b2;
79 jgs 150 }
80     }
81     return;
82     }

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26