Actual source code: taomat_petsc.c

  1: #include "tao_general.h"

  3: #ifdef TAO_USE_PETSC

  5: #include "taomat_petsc.h"
  6: #include "../vector/taovec_petsc.h"
  7: #include "../indexset/taois_petsc.h"

  9: int MatD_Fischer(Mat, Vec, Vec, Vec, Vec, Vec, Vec, Vec, Vec);
 10: int MatD_SFischer(Mat, Vec, Vec, Vec, Vec, double, Vec, Vec, Vec, Vec, Vec);
 11: int MatSMFResetRowColumn(Mat,IS,IS);

 15: /*@C
 16:    TaoWrapPetscMat - Creates a new TaoMat object using a PETSc matrix.

 18:    Input Parameter:
 19: +  M -  a PETSc matrix
 20: -  MM -  the address of a pointer to a TaoMatPetsc

 22:    Output Parameter:
 23: .  MM - the address of a pointer to new TaoMat

 25:    Note:  
 26:    A TaoMatPetsc is an object with the methods of an abstract
 27:    TaoMat object.  A TaoMatPetsc contains an implementation of the TaoMat
 28:    methods.  Routines using these vectors should declare a pointer to 
 29:    a TaoMat, assign this pointer to the address of a TaoMat object, 
 30:    use the pointer to invoke methods on the object, and use this pointer
 31:    as an argument when calling other routines.  This usage is different
 32:    from the usage of a PETSc Mat.  In PETSc, applications will typically
 33:    declare a Mat, and pass it as an argument into routines.  That is,
 34:    applications will typically declare a pointer to a TaoMat and use the
 35:    pointer, or declare a Mat and use it directly.

 37:    Note:
 38:    The user is repsonsible for destroying the Mat M, in addition to
 39:    to TaoMatPetsc vector MM.  The Mat can be destroyed immediately
 40:    after this routine.

 42:   Level: developer

 44: .seealso TaoMatGetPetscMat(), TaoMatDestroy()
 45: @*/
 46: int TaoWrapPetscMat( Mat M, TaoMatPetsc* *MM){
 48:   if (MM){ *MM = new  TaoMatPetsc(M);}
 49:   return(0);
 50: }

 54: /* @C
 55:    TaoMatPetsc - Creates a new TaoMat object using a PETSc matrix.

 57:    Input Parameter:
 58: +  MM -  a PETSc matrix

 60:    Level: advanced

 62:    Note:
 63:    The method TaoMatPetsc::SetMatrix(Mat M, Mat MPre, MatStructure) replaces the matrix MM with M.
 64:    The method TaoMatPetsc::GetMatrix() returns a pointer to the matrix MM, its preconditioner, and KSP flag

 66: @ */
 67: TaoMatPetsc::TaoMatPetsc(Mat MM) : TaoMat()
 68: {
 69:   int info;
 70:   pm=0; pm_pre=0;
 71:   preflag=DIFFERENT_NONZERO_PATTERN;
 72:   MPI_Comm comm;
 73:   if (MM) { 
 74:     PetscObjectGetComm((PetscObject)MM,&comm);
 75:     pmatviewer=PETSC_VIEWER_STDOUT_(comm);
 76:     SetMatrix(MM,MM,preflag);  
 77:     info = PetscInfo(MM,"Wrap a PETSc Mat to create a TaoMat .\n");
 78:   }
 79:   return;
 80: }

 84: TaoMatPetsc::~TaoMatPetsc()
 85: {
 86:   int info;
 87:   if (pm) {
 88:     info = PetscInfo(pm,"TAO: Destroy a PETSc Mat within a TaoMat .\n");
 89:     PetscObjectDereference((PetscObject)pm);
 90:     PetscObjectDereference((PetscObject)pm_pre);
 91:   }
 92:   pm=0; pm_pre=0;
 93:   return;
 94: }

 98: /* @C
 99:    GetMatrix - Gets the underlying PETSc matrix information

101:    Output Parameters:
102: +  M -  a PETSc matrix
103: .  MPre -  a PETSc preconditioning matrix (use PETSC_NULL for no change)
104: -  flag -  flag associated with KSPSetOperators (PETSC_NULL for no change).

106:    Level: advanced

108: .seealso TaoMatPetsc::SetMatrix()
109:  
110: @ */
111: int TaoMatPetsc::GetMatrix(Mat *M, Mat *MPre, MatStructure *flag){
113:   if (M) *M = pm;
114:   if (flag) *flag = preflag;
115:   if (MPre) *MPre = pm_pre;
116:   return(0);  
117: }


122: /* @C
123:    SetMatrix - Changes the underlying PETSc matrix structrure

125:    Input Parameters:
126: +  M -  a PETSc matrix
127: .  MPre -  a PETSc preconditioning matrix (use PETSC_NULL for no change)
128: -  flag -  flag associated with KSPSetOperators (PETSC_NULL for no change).

130:    Level: advanced

132:    Note:
133:    The method TaoMatPetsc::GetMatrix() returns a pointer to the matrix M

135: @ */
136: int TaoMatPetsc::SetMatrix(Mat M, Mat MPre, MatStructure flag){
137:   int info;
138:   int r1,r2,c1,c2,R1,R2,C1,C2;
140:   if (M || MPre){
143:     if (M!=MPre){
144:       info = MatGetSize(M,&R1,&C1);CHKERRQ(info);
145:       info = MatGetSize(MPre,&R2, &C2);CHKERRQ(info);
146:       info = MatGetLocalSize(M,&r1,&c1);CHKERRQ(info);
147:       info = MatGetLocalSize(MPre,&r2,&c2);CHKERRQ(info);
148:       if (R1!=R2 || C1!=C2 || r1!=r2 || c1!=c2) {
149:         SETERRQ(1,"TAO Error: Preconditioning Matrix does not match the Solution matrix");
150:       }
151:     }
152:     
153:     info = PetscObjectReference((PetscObject)M);CHKERRQ(info);
154:     info = PetscObjectReference((PetscObject)MPre);CHKERRQ(info); 
155:   }

157:   if (pm&&pm_pre){
158:     info = PetscObjectDereference((PetscObject)pm);CHKERRQ(info);
159:     info = PetscObjectDereference((PetscObject)pm_pre);CHKERRQ(info); 
160:   }

162:   pm=M;
163:   pm_pre=MPre; 
164:   preflag=flag;
165:   return(0);  
166: }



172: int TaoMatPetsc::Compatible(TaoVec* xx, TaoVec* yy, TaoTruth *flag){
173:   int info;
174:   int pmCol,pmRow,pmrow,pmcol,xN,yN,yn;
175:   Vec x,y;

178:   *flag=TAO_TRUE;
179:   if (xx==0 || yy==0){
180:     *flag=TAO_FALSE;
181:     return(0);
182:   }
183:   x=((TaoVecPetsc*)xx)->GetVec();
184:   y=((TaoVecPetsc*)yy)->GetVec();
185:   info = MatGetSize(pm,&pmRow,&pmCol); CHKERRQ(info);
186:   info = MatGetLocalSize(pm,&pmrow,&pmcol); CHKERRQ(info);
187:   info = VecGetSize(x,&xN); CHKERRQ(info);
188:   info = VecGetSize(y,&yN); CHKERRQ(info);
189:   info = VecGetLocalSize(y,&yn); CHKERRQ(info);
190:   if ((pmCol != xN) || (pmRow != yN) || (pmrow != yn)) {
191:     *flag=TAO_FALSE;
192:     return(0);
193:   }
194:   
198:   
199:   return(0);
200: }

204: int TaoMatPetsc::Clone(TaoMat* *ntm){
205:   int info;
206:   TaoMatPetsc *nptm;
207:   Mat M,MPre,npm,nmpre;
208:   MatStructure flag;

211:   info = GetMatrix(&M,&MPre,&flag); CHKERRQ(info);
212:   info = MatDuplicate(M,MAT_COPY_VALUES,&npm); CHKERRQ(info);
213:   info = MatAssemblyBegin(npm,MAT_FINAL_ASSEMBLY); CHKERRQ(info);
214:   info = MatAssemblyEnd(npm,MAT_FINAL_ASSEMBLY); CHKERRQ(info);
215:   if (M!=MPre){
216:     info = MatDuplicate(MPre,MAT_COPY_VALUES,&nmpre); CHKERRQ(info);
217:     info = MatAssemblyBegin(nmpre,MAT_FINAL_ASSEMBLY); CHKERRQ(info);
218:     info = MatAssemblyEnd(nmpre,MAT_FINAL_ASSEMBLY); CHKERRQ(info);
219:   }
220:   info = TaoWrapPetscMat(npm, &nptm);   CHKERRQ(info);
221:   *ntm = nptm;
222:   info = nptm->SetMatrix(npm,nmpre,flag); CHKERRQ(info);
223:   info = MatDestroy(npm); CHKERRQ(info);
224:   if (M!=MPre){
225:     info = MatDestroy(nmpre); CHKERRQ(info);
226:   }
227:   return(0);
228: }

232: int TaoMatPetsc::CopyFrom(TaoMat* tm){
233:   int info;
234:   TaoMatPetsc* tmm=(TaoMatPetsc*)tm;
235:   Mat M,MPre;
236:   MatStructure flag;

239:   info = tmm->GetMatrix(&M,&MPre,&flag);CHKERRQ(info);
240:   info = MatCopy(M,pm,SAME_NONZERO_PATTERN); 
241:   CHKERRQ(info);
242:   if (pm!=pm_pre){
243:     info = MatCopy(MPre,pm_pre,SAME_NONZERO_PATTERN); CHKERRQ(info);
244:   }
245:   return(0);
246: }

250: int TaoMatPetsc::GetDimensions( int *m, int *n ){
251:   int info;
253:   info=MatGetSize(pm,m,n); CHKERRQ(info);
254:   return(0);  
255: }

259: int TaoMatPetsc::Multiply(TaoVec *tv, TaoVec *tw)
260: {
261:   TaoVecPetsc *pv = dynamic_cast <TaoVecPetsc *> (tv);
262:   TaoVecPetsc *pw = dynamic_cast <TaoVecPetsc *> (tw);

264:   int info;
266:   info=MatMult(pm,pv->GetVec(),pw->GetVec()); CHKERRQ(info);
267:   return(0);
268: }

272: int TaoMatPetsc::MultiplyAdd(TaoVec* tv,TaoVec* tw,TaoVec* ty)
273: {
274:   TaoVecPetsc *pv = dynamic_cast <TaoVecPetsc *> (tv);
275:   TaoVecPetsc *pw = dynamic_cast <TaoVecPetsc *> (tw);
276:   TaoVecPetsc *py = dynamic_cast <TaoVecPetsc *> (ty);
277:   int info;
279:   info=MatMultAdd(pm,pv->GetVec(),pw->GetVec(),py->GetVec()); CHKERRQ(info);
280:   return(0);
281: }

285: int TaoMatPetsc::MultiplyTranspose(TaoVec* tv,TaoVec* tw)
286: {
287:   TaoVecPetsc *pv = dynamic_cast <TaoVecPetsc *> (tv);
288:   TaoVecPetsc *pw = dynamic_cast <TaoVecPetsc *> (tw);
289:   int info;
291:   info=MatMultTranspose(pm,pv->GetVec(),pw->GetVec()); CHKERRQ(info);
292:   return(0);
293: }

297: int TaoMatPetsc::MultiplyTransposeAdd(TaoVec* tv,TaoVec* tw,TaoVec* ty)
298: {
299:   TaoVecPetsc *pv = dynamic_cast <TaoVecPetsc *> (tv);
300:   TaoVecPetsc *pw = dynamic_cast <TaoVecPetsc *> (tw);
301:   TaoVecPetsc *py = dynamic_cast <TaoVecPetsc *> (ty);
302:   int info;
304:   info=MatMultTransposeAdd(pm,pv->GetVec(),pw->GetVec(),py->GetVec()); CHKERRQ(info);
305:   return(0);
306: }

310: int TaoMatPetsc::SetDiagonal(TaoVec* tv)
311: {
312:   TaoVecPetsc *pv = dynamic_cast <TaoVecPetsc *> (tv);
313:   int info;
315:   info = MatDiagonalSet(pm,pv->GetVec(),INSERT_VALUES); CHKERRQ(info);
316:   if (pm != pm_pre){
317:     info = MatDiagonalSet(pm_pre,pv->GetVec(),INSERT_VALUES); CHKERRQ(info);
318:   }
319:   return(0);
320: }

324: int TaoMatPetsc::AddDiagonal(TaoVec* tv)
325: {
326:   TaoVecPetsc *pv = dynamic_cast <TaoVecPetsc *> (tv);
327:   int info;
329:   info = MatDiagonalSet(pm,pv->GetVec(),ADD_VALUES); CHKERRQ(info);
330:   if (pm != pm_pre){
331:     info = MatDiagonalSet(pm_pre,pv->GetVec(),ADD_VALUES); CHKERRQ(info);
332:   }
333:   return(0);
334: }

338: int TaoMatPetsc::GetDiagonal(TaoVec* tv)
339: {
340:   TaoVecPetsc *pv = dynamic_cast <TaoVecPetsc *> (tv);
341:   int info;
343:   info = MatGetDiagonal(pm,pv->GetVec()); CHKERRQ(info);
344:   return(0);
345: }

349: int TaoMatPetsc::ShiftDiagonal(double c)
350: {
351:   int info;
352:   PetscScalar cc=c;
354:   info = MatShift(pm,cc); CHKERRQ(info);
355:   if (pm != pm_pre) {
356:     info = MatShift(pm_pre,cc); CHKERRQ(info);
357:   }
358:   return(0);
359: }

363: /*@C
364:    SetPetscViewer

366:    Input Parameter:
367: .  viewer - a viewer object to be used with View() and MatView()

369:    Level: advanced

371: @*/
372: int TaoMatPetsc::SetPetscViewer(PetscViewer viewer){
374:   pmatviewer= viewer;
375:   return(0);  
376: }


381: int TaoMatPetsc::View()
382: {
383:   int info;
384:   PetscTruth flg=PETSC_FALSE;
386:   info = PetscOptionsHasName(PETSC_NULL,"-tao_mat_draw",&flg); CHKERRQ(info);
387:   if (flg){
388:     info = MatView(pm,PETSC_VIEWER_DRAW_WORLD); CHKERRQ(info);
389:   } else {
390:     info = MatView(pm,pmatviewer); CHKERRQ(info);
391:   }
392:   return(0);
393: }


398: int TaoMatPetsc::RowScale(TaoVec* tv)
399: {
400:   TaoVecPetsc *pv = dynamic_cast <TaoVecPetsc *> (tv);
401:   int info;

404:   info = MatDiagonalScale(pm, pv->GetVec(), PETSC_NULL); CHKERRQ(info);
405:   if (pm!=pm_pre){
406:     info = MatDiagonalScale(pm_pre,pv->GetVec(),PETSC_NULL); CHKERRQ(info);
407:   }
408:   return(0);
409: }

413: int TaoMatPetsc::ColScale(TaoVec* tv){
414:   Vec vv = ((TaoVecPetsc *)tv)->GetVec();
415:   int info;

418:   info = MatDiagonalScale(pm, PETSC_NULL, vv); CHKERRQ(info);
419:   if (pm!=pm_pre){
420:     info = MatDiagonalScale(pm_pre,PETSC_NULL,vv); CHKERRQ(info);
421:   }
422:   return(0);
423: }

427: int TaoMatPetsc::CreateReducedMatrix(TaoIndexSet* S1,TaoIndexSet* S2,TaoMat* *MM){
428:   TaoIndexSetPetsc *prows = dynamic_cast <TaoIndexSetPetsc *> (S1);
429:   TaoIndexSetPetsc *pcols = dynamic_cast <TaoIndexSetPetsc *> (S2);
430:   int info;
431:   TaoMatPetsc *M;
432:   Mat B=0,BPre,A,APre;
433:   TaoPetscISType type;
434:   MatStructure flag;
435:   PetscTruth assembled,flg;

438:   info = prows->GetReducedType(&type); CHKERRQ(info);
439:   if (type==TaoMaskFullSpace){
440:     info = GetMatrix(&B,&BPre,&flag); CHKERRQ(info);
441:     PetscOptionsHasName(0,"-different_submatrix",&flg);
442:     if (flg==PETSC_TRUE){
443:       info = TaoWrapPetscMat(B, &M); CHKERRQ(info);
444:       info = MatAssembled(B,&assembled); CHKERRQ(info);
445:       if (assembled==PETSC_TRUE){
446:         info = MatDuplicate(B,MAT_DO_NOT_COPY_VALUES,&A); CHKERRQ(info);
447:         info = M->SetMatrix(A,A,flag); CHKERRQ(info);
448:         info = MatDestroy(A); CHKERRQ(info);
449:       }
450:       info = MatAssembled(BPre,&assembled); CHKERRQ(info);
451:       if (B != BPre && assembled==PETSC_TRUE){
452:         info = MatDuplicate(BPre,MAT_DO_NOT_COPY_VALUES,&APre); CHKERRQ(info);
453:         info = M->SetMatrix(A,APre,flag); CHKERRQ(info);
454:         info = MatDestroy(APre); CHKERRQ(info);
455:       }
456:     } else {
457:       info = GetMatrix(&B,&BPre,&flag); CHKERRQ(info);
458:       info = TaoWrapPetscMat(B, &M); CHKERRQ(info);
459:       info = M->SetMatrix(B,BPre,flag); CHKERRQ(info);
460:     }    
461:   } else if (type==TaoMatrixFree){
462:     info = GetMatrix(&B,&BPre,&flag); CHKERRQ(info);
463:     info = MatCreateSubMatrixFree(B,prows->GetIS(), pcols->GetIS(), &A); CHKERRQ(info);
464:     info = TaoWrapPetscMat(A, &M); CHKERRQ(info);
465:     info = M->SetMatrix(A,BPre,flag); CHKERRQ(info);
466:     info = PetscObjectDereference((PetscObject)A);CHKERRQ(info);
467:   } else {
468:     info = GetMatrix(&B,&BPre,&flag); CHKERRQ(info);
469:     info = TaoWrapPetscMat(B, &M); CHKERRQ(info);
470:     info = M->SetMatrix(B,BPre,flag); CHKERRQ(info);
471:   }
472:   *MM=M;

474:   return(0);
475: }

479: int TaoMatPetsc::SetReducedMatrix(TaoMat* M,TaoIndexSet* S1,TaoIndexSet* S2)
480: {
481:   TaoIndexSetPetsc *prows = dynamic_cast <TaoIndexSetPetsc *> (S1);
482:   TaoIndexSetPetsc *pcols = dynamic_cast <TaoIndexSetPetsc *> (S2);

484:   int info;
485:   TaoPetscISType type;
486:   TaoMatPetsc *MM=((TaoMatPetsc *)M);
487:   Mat A,Apre,B,BPre;
488:   MatStructure flag;

491:   info=GetMatrix(&B,&BPre,0); CHKERRQ(info);
492:   info=MM->GetMatrix(&A,&Apre,0); CHKERRQ(info);
493:   info = prows->GetReducedType(&type); CHKERRQ(info);
494:   if (type==TaoMaskFullSpace){
495:     Vec DDiag,VRow,VCol;
496:     info = GetMatrix(&B,&BPre,&flag); CHKERRQ(info);
497:     info = MM->GetMatrix(&A,&Apre,&flag); CHKERRQ(info);
498:     if (!A){
499:       info = MatDuplicate(B,MAT_COPY_VALUES,&A); CHKERRQ(info);    
500:       if (B!=BPre){info = MatDuplicate(BPre,MAT_COPY_VALUES,&Apre); CHKERRQ(info);}
501:       info = MM->SetMatrix(A,Apre,flag); CHKERRQ(info);
502:       info=MatDestroy(A); CHKERRQ(info);
503:       if (B!=BPre){info=MatDestroy(Apre); CHKERRQ(info);}
504:     }
505:     if (A!=B){ info=MatCopy(A,B,SAME_NONZERO_PATTERN); CHKERRQ(info); }
506:     if (B!=BPre && Apre!=BPre){ info=MatCopy(Apre,BPre,SAME_NONZERO_PATTERN); CHKERRQ(info); }
507:     info=prows->GetMask(&VRow); CHKERRQ(info);
508:     info=pcols->GetMask(&VCol); CHKERRQ(info);
509:     info=VecDuplicate(VRow,&DDiag); CHKERRQ(info);
510:     info=MatGetDiagonal(A,DDiag); CHKERRQ(info);
511:     info=MatDiagonalScale(B,VRow,VCol); CHKERRQ(info);
512:     info=MatDiagonalSet(B,DDiag,INSERT_VALUES); CHKERRQ(info);
513:     if (B!=BPre){
514:       info=MatGetDiagonal(Apre,DDiag); CHKERRQ(info);
515:       info=MatDiagonalScale(BPre,VRow,VCol); CHKERRQ(info);
516:       info=MatDiagonalSet(BPre,DDiag,INSERT_VALUES); CHKERRQ(info);
517:     }
518:     info=VecDestroy(DDiag); CHKERRQ(info);
519:   } else if (type==TaoMatrixFree){
520:     info = MM->GetMatrix(&A,&Apre,&flag); CHKERRQ(info);
521:     info = MatCreateSubMatrixFree(A,prows->GetIS(),pcols->GetIS(),&B); CHKERRQ(info);
522:     info = SetMatrix(B,Apre,flag); CHKERRQ(info);
523:     info = PetscObjectDereference((PetscObject)B);CHKERRQ(info);
524:   } else {
525:     IS LocalRows, LocalCols, AllCols;
526:     int nn;
527:     //    info = GetMatrix(&B,&BPre,&flag); CHKERRQ(info);
528:     info = MM->GetMatrix(&A,&Apre,&flag); CHKERRQ(info);
529:     info = prows->RedistributeIS(&LocalRows); CHKERRQ(info);
530:     info = pcols->RedistributeIS(&LocalCols); CHKERRQ(info);
531:     info = pcols->GetWholeIS(&AllCols); CHKERRQ(info);
532:     info = ISGetLocalSize(LocalCols,&nn); CHKERRQ(info);
533:     info = MatGetSubMatrix(A,LocalRows,AllCols,nn,MAT_INITIAL_MATRIX,&B);
534:     CHKERRQ(info);
535:     BPre=B;

537:     if (A != Apre){      
538:       info = MatGetSubMatrix(Apre,LocalRows,AllCols,nn,MAT_INITIAL_MATRIX,&BPre);
539:       CHKERRQ(info);
540:     } else {
541:       info = PetscObjectReference((PetscObject)BPre);CHKERRQ(info);
542:     }
543:     
544:     info = SetMatrix(B,BPre,flag); CHKERRQ(info);
545:     info = PetscObjectDereference((PetscObject)B);CHKERRQ(info);
546:     info = PetscObjectDereference((PetscObject)BPre);CHKERRQ(info);
547:   }

549:   return(0);
550: }

554: int TaoMatPetsc::D_Fischer(TaoVec* tx, TaoVec* tf, TaoVec* tl, TaoVec* tu,
555:                            TaoVec* tt1, TaoVec* tt2, TaoVec* tda, TaoVec* tdb) 
556: {
557:   Vec x = ((TaoVecPetsc *)tx)->GetVec();
558:   Vec f = ((TaoVecPetsc *)tf)->GetVec();
559:   Vec l = ((TaoVecPetsc *)tl)->GetVec();
560:   Vec u = ((TaoVecPetsc *)tu)->GetVec();
561:   Vec da = ((TaoVecPetsc *)tda)->GetVec();
562:   Vec db = ((TaoVecPetsc *)tdb)->GetVec();
563:   Vec t1 = ((TaoVecPetsc *)tt1)->GetVec();
564:   Vec t2 = ((TaoVecPetsc *)tt2)->GetVec();
565:   int info;

568:   info = MatD_Fischer(pm, x, f, l, u, t1, t2, da, db); CHKERRQ(info);
569:   return(0);
570: }

574: int TaoMatPetsc::D_SFischer(TaoVec* tx, TaoVec* tf, TaoVec* tl, TaoVec* tu,
575:                             double mu,
576:                             TaoVec* tt1, TaoVec* tt2,
577:                             TaoVec* tda, TaoVec* tdb, TaoVec* tdm)
578: {
579:   Vec x = ((TaoVecPetsc *)tx)->GetVec();
580:   Vec f = ((TaoVecPetsc *)tf)->GetVec();
581:   Vec l = ((TaoVecPetsc *)tl)->GetVec();
582:   Vec u = ((TaoVecPetsc *)tu)->GetVec();
583:   Vec t1 = ((TaoVecPetsc *)tt1)->GetVec();
584:   Vec t2 = ((TaoVecPetsc *)tt2)->GetVec();
585:   Vec da = ((TaoVecPetsc *)tda)->GetVec();
586:   Vec db = ((TaoVecPetsc *)tdb)->GetVec();
587:   Vec dm = ((TaoVecPetsc *)tdm)->GetVec();
588:   int info;

591:   if ((mu >= -TAO_EPSILON) && (mu <= TAO_EPSILON)) {
592:     tdm->SetToZero();
593:     D_Fischer(tx, tf, tl, tu, tt1, tt2, tda, tdb);
594:   }
595:   else {
596:     info = MatD_SFischer(pm, x, f, l, u, mu, t1, t2, da, db, dm); CHKERRQ(info);
597:   }
598:   return(0);
599: }

603: int TaoMatPetsc::Norm1(double *nm){
604:   int info;
605:   PetscReal nnmm;
607:   info = MatNorm(pm,NORM_1,&nnmm);CHKERRQ(info);
608:   *nm=nnmm;
609:   return(0);
610: }

612: #endif