Actual source code: taomat_ga.c
1: //File: taomat_ga.c
3: /**************************************************************
5: Author: Limin Zhang, Ph.D.
6: Mathematics Department
7: Columbia Basin College
8: Pasco, WA 99301
9: Limin.Zhang@cbc2.org
11: Mentor: Jarek Naplocha, Ph.D.
12: Environmental Molecular Science Laboratory
13: Pacific Northwest National Laboratory
14: Richland, WA 99352
16: Date: 7/11/2002
18: Purpose:
19: to design and implement some interfaces between TAO and
20: global arrays.
22: Revise History:
24: 8/8/02
25: To clean Petsc function calls and marcos.
26: **************************************************************/
29: //#include "tao_general.h"
31: #include "taomat_ga.h"
32: #include "taovec_ga.h"
34: /*
35: int MatD_Fischer(Mat, Vec, Vec, Vec, Vec, Vec, Vec, Vec, Vec);
36: int MatCreateADA(Mat,Vec,Vec,Mat*);
37: int MatSMFResetRowColumn(Mat,IS,IS);
38: */
40: /*@C
41: TaoWrapGaMat - Creates a new TaoMat object using a GA matrix.
43: Input Parameter:
44: . M - a GA matrix
46: Output Parameter:
47: . MM - new TaoMat
49: Level: advanced
51: .seealso TaoMatGetGaMat(), TaoMatDestroy()
52: @*/
53: int TaoWrapGaMat (GAMat M, TaoMatGa ** MM)
54: {
55: TaoFunctionBegin;
56: *MM = new TaoMatGa (M);
57: TaoFunctionReturn (0);
58: }
60: /*@C
61: TaoMatGetGaMat - If the TaoMat is of the TaoMatGa type, this routine
62: gets the underlying GA matrix.
64: Input Parameter:
65: . MM - the TaoMat
67: Output Parameter:
68: . M - the GA mat
70: Note:
71: The function TaoMatGa::GetMat() will also return the Mat M.
73: Level: advanced
74: @*/
75: int TaoMatGetGaMat (TaoMat * MM, GAMat * M)
76: {
77: TaoFunctionBegin;
78: if (MM)
79: {
80: *M = ((TaoMatGa *) MM)->GetMat ();
81: }
82: else
83: {
84: *M = 0;
85: }
86: TaoFunctionReturn (0);
87: }
89: TaoMatGa::TaoMatGa (GAMat MM):TaoMat ()
90: {
91: this->pm = MM;
92: this->MatObject = (void *) (&MM);
93: this->pm_pre = MM;
94: return;
95: }
97: int TaoMatGa::Compatible (TaoVec * xx, TaoVec * yy, TaoTruth *flag)
98: {
99: int
100: nx,
101: ny,
102: m,
103: n;
104: GAVec
105: x,
106: y;
107: if (xx == 0 || yy == 0)
108: { *flag=TAO_FALSE; return 0; }
110: x = ((TaoVecGa *) xx)->GetVec ();
111: y = ((TaoVecGa *) yy)->GetVec ();
112: if (x == 0 || y == 0 || this->pm == 0)
113: { *flag=TAO_FALSE; return 0; }
115: int
116: xtype,
117: ytype,
118: type,
119: ndim,
120: dims[2];
122: NGA_Inquire (x, &xtype, &ndim, &nx);
123: if (ndim != 1)
124: GA_Error ("TaoMatGa::Compatible: Wrong dimension", ndim);
126: NGA_Inquire (y, &ytype, &ndim, &ny);
127: if (ndim != 1)
128: GA_Error ("TaoMatGa::Compatible:Wrong dimension", ndim);
130: NGA_Inquire (this->pm, &type, &ndim, dims);
131: if (ndim != 2)
132: GA_Error ("TaoMatGa::Compatible:Wrong dimension", ndim);
133: n = dims[0];
134: m = dims[1];
136: if (xtype != type || ytype != type)
137: GA_Error ("TaoMatGa::Compatible:Wrong data type", type);
139: if (n != nx || m != ny)
140: { *flag=TAO_FALSE; return 0; }
142: *flag=TAO_TRUE;
143: return 0;
145: }
147: int
148: TaoMatGa::Clone (TaoMat ** ntm)
149: {
150: TaoMatGa *nptm;
151: GAMat npm;
152: TaoFunctionBegin;
153: npm = GA_Duplicate (this->pm, "Cloned");
154: if (!npm)
155: GA_Error ("TaoMatGa::Clone:duplicate failed: ", npm);
156: TaoWrapGaMat (npm, &nptm);
157: *ntm = nptm;
158: TaoFunctionReturn (0);
159: }
161: int
162: TaoMatGa::CopyFrom (TaoMat * tm)
163: {
164: TaoFunctionBegin;
165: #if 0 //try a different way to get GAMat
166: TaoMatGa *temp2 = dynamic_cast < TaoMatGa * >(tm);
167: GAMat pm2 = temp2->GetMat ();
168: #endif
169: GAMat pm2 = ((TaoMatGa *)tm)->GetMat ();
170: GA_Copy (pm2, this->pm);
171: TaoFunctionReturn (0);
172: }
174: int
175: TaoMatGa::GetDimensions (int *m, int *n)
176: {
177: TaoFunctionBegin;
178: int type, ndim, dims[2];
179: NGA_Inquire (this->pm, &type, &ndim, dims);
180: if (ndim != 2)
181: GA_Error ("TaoMatGa::GetDimension: wrong dimension", ndim);
182: else
183: {
184: *n = dims[0];
185: *m = dims[1];
186: }
187: TaoFunctionReturn (0);
188: }
190: int
191: TaoMatGa::Multiply (TaoVec * tv, TaoVec * tw)
192: {
193: double alpha = 1.0;
194: double beta = 0.0;
195: TaoTruth flag;
196: TaoFunctionBegin;
197: GAVec vv = ((TaoVecGa *) tv)->GetVec ();
198: GAVec ww = ((TaoVecGa *) tw)->GetVec ();
199: int m, n = 1, k; //n = 1 due to vecotor ww that is k by 1 array
200: this->GetDimensions (&m, &k);
201: this->Compatible (tv, tw,&flag);
202: if (flag)
203: GA_Dgemm ('N', 'N', m, n, k, alpha, this->pm, vv, beta, ww);
204: else
205: TaoFunctionReturn (1);
206: TaoFunctionReturn (0);
207: }
209: int
210: TaoMatGa::MultiplyAdd (TaoVec * tv, TaoVec * tw, TaoVec * ty)
211: {
212: TaoTruth flag1,flag2;
213: TaoFunctionBegin;
214: GAVec xx = ((TaoVecGa *) tv)->GetVec ();
215: GAVec ww = ((TaoVecGa *) tw)->GetVec ();
216: GAVec yy = ((TaoVecGa *) ty)->GetVec ();
217: GAVec temp_ww = GA_Duplicate (ww, "temp_WW");
218: if (!temp_ww)
219: GA_Error ("TaoMatGa::MultiplyAdd:Duplicate:", temp_ww);
220: GA_Copy (ww, temp_ww);
221: double alpha = 1.0;
222: double beta = 1.0;
223: int m, n = 1, k; //n = 1 due to vecotor ww that is k by 1 array
224: this->GetDimensions (&m, &k);
225: this->Compatible (tv, ty, &flag1);
226: this->Compatible (tv, tw, &flag2);
227: if (flag1 && flag2)
228: {
229: GA_Dgemm ('N', 'N', m, n, k, alpha, this->pm, xx, beta, temp_ww);
230: //copy from temp_ww to yy. Doing this way, the data in ww is not dirty!!!!!
231: GA_Copy (temp_ww, yy);
232: //free the temporary memory
233: GA_Destroy(temp_ww);
234: }
235: else
236: TaoFunctionReturn (1);
237: TaoFunctionReturn (0);
238: }
240: int
241: TaoMatGa::MultiplyTranspose (TaoVec * tv, TaoVec * tw)
242: {
243: double alpha = 1.0;
244: double beta = 0.0;
245: int m, n = 1, k; /*n = 1 due to vecotor ww that is k */
246: TaoFunctionBegin;
247: GAVec vv = ((TaoVecGa *) tv)->GetVec ();
248: GAVec ww = ((TaoVecGa *) tw)->GetVec ();
249: this->GetDimensions (&m, &k);
250: GA_Dgemm ('T', 'N', k, n, m, alpha, this->pm, vv, beta, ww);
251: TaoFunctionReturn (0);
252: }
254: int
255: TaoMatGa::MultiplyTransposeAdd (TaoVec * tv, TaoVec * tw, TaoVec * ty)
256: {
257: GAVec xx = ((TaoVecGa *) tv)->GetVec ();
258: GAVec ww = ((TaoVecGa *) tw)->GetVec ();
259: GAVec yy = ((TaoVecGa *) ty)->GetVec ();
260: TaoFunctionBegin;
261: TaoFunctionBegin;
262: GAVec temp_ww = GA_Duplicate (ww, "temp_WW");
263: if (!temp_ww)
264: GA_Error ("TaoMatGa::MultiplyAdd:Dupli cate:", temp_ww);
265: GA_Copy (ww, temp_ww);
266: double alpha = 1.0;
267: double beta = 1.0;
268: int m, n = 1, k;
269: //n = 1 due to vecotor ww that is k by 1 array
270: this->GetDimensions (&m, &k);
271: GA_Dgemm ('T', 'N', k, n, m, alpha, this->pm, xx, beta, temp_ww);
272: //copy from temp_ww to yy. Doing this way, the data in ww is not dirty!!!!!
273: GA_Copy (temp_ww, yy);
274: //free the memory held by the temporary GA array temp_ww
275: GA_Destroy(temp_ww);
276: TaoFunctionReturn (0);
277: }
279: int
280: TaoMatGa::SetDiagonal (TaoVec * tv)
281: {
282: GAVec vv = ((TaoVecGa *) tv)->GetVec ();
283: TaoFunctionBegin;
284: GA_Set_diagonal (this->pm, vv);
285: TaoFunctionReturn (0);
286: }
288: int
289: TaoMatGa::AddDiagonal (TaoVec * tv)
290: {
291: GAVec vv = ((TaoVecGa *) tv)->GetVec ();
292: TaoFunctionBegin;
293: GA_Add_diagonal (this->pm, vv);
294: TaoFunctionReturn (0);
295: }
297: int
298: TaoMatGa::GetDiagonal (TaoVec * tv)
299: {
300: GAVec vv = ((TaoVecGa *) tv)->GetVec ();
301: TaoFunctionBegin;
302: GA_Get_diag (this->pm, vv);
303: TaoFunctionReturn (0);
304: }
306: int
307: TaoMatGa::ShiftDiagonal (TaoScalar c)
308: {
309: TaoFunctionBegin;
310: GA_Shift_diagonal (this->pm, &c);
311: TaoFunctionReturn (0);
312: }
314: int
315: TaoMatGa::View ()
316: {
317: TaoFunctionBegin;
318: GA_Print (this->pm);
319: TaoFunctionReturn (0);
320: }
323: int
324: TaoMatGa::Presolve ()
325: {
326: TaoFunctionBegin;
327: TaoFunctionReturn (0);
328: }
330: int
331: TaoMatGa::Solve (TaoVec * tv, TaoVec * tw, TaoTruth * tt)
332: {
333: TaoFunctionBegin;
334: GA_Error ("not yet implemented.", 0);
335: TaoFunctionReturn (1);
336: }
338: int
339: TaoMatGa::RowScale (TaoVec * tv)
340: {
341: GAVec vv = ((TaoVecGa *) tv)->GetVec ();
342: TaoFunctionBegin;
343: GA_Scale_rows (this->pm, vv);
344: TaoFunctionReturn (0);
345: }
347: int
348: TaoMatGa::ColScale (TaoVec * tv)
349: {
350: GAVec vv = ((TaoVecGa *) tv)->GetVec ();
351: TaoFunctionBegin;
352: GA_Scale_cols (this->pm, vv);
353: TaoFunctionReturn (0);
354: }
356: int
357: TaoMatGa::CreateReducedMatrix (TaoIndexSet * S1,
358: TaoIndexSet * S2, TaoMat ** MM)
359: {
360: int info;
361: TaoMatGa *M;
362: GAMat B = 0;
363: TaoFunctionBegin;
364: info = TaoWrapGaMat (B, &M);
365: CHKERRQ (info);
366: *MM = M;
367: TaoFunctionReturn (0);
368: }
370: int
371: TaoMatGa::SetReducedMatrix (TaoMat * M, TaoIndexSet * S1, TaoIndexSet * S2)
372: {
373: TaoFunctionBegin;
374: GA_Error ("Not yet implemented.", 0);
375: TaoFunctionReturn (1);
376: }
380: int
381: TaoMatGa::D_Fischer (TaoVec * tx, TaoVec * tf, TaoVec * tl,
382: TaoVec * tu, TaoVec * tda,
383: TaoVec * tdb, TaoVec * tt1, TaoVec * tt2)
384: {
385: TaoFunctionBegin;
386: GA_Error ("Not yet implemented.", 0);
387: TaoFunctionReturn (1);
388: }
390: int
391: TaoMatGa::Norm1 (double *nm)
392: {
393: TaoFunctionBegin;
394: GA_Norm1 (this->pm, nm);
395: TaoFunctionReturn (0);
396: }
399: int
400: TaoMatGa::Fill (double c)
401: {
403: TaoFunctionBegin;
404: GA_Fill (this->pm, &c);
405: TaoFunctionReturn (0);
406: }
409: int
410: TaoMatGa::Zero ()
411: {
413: TaoFunctionBegin;
414: GA_Zero (this->pm);
415: TaoFunctionReturn (0);
417: }