Actual source code: daapp_solve.c
1: #include "taodaapplication.h" /*I "taodaapplication.h" I*/
2: //#include "taoapp.h"
3: #include "tao.h"
4: #include "petsc.h"
5: #include "src/petsctao/application/petscapp/tao_app_impl.h"
6: #include "daapp_impl.h"
8: static int Tao_DA_Solve=0;
10: int TaoAppDAApp(TAO_APPLICATION, DA_APPLICATION *);
14: /*@
15: TaoDAAppSolve - Solve the PETSC DA application.
17: Input Parameters:
18: . daapplication - the TAO DAApplication structure
20: Level: beginner
22: .keywords: Application, DA, Solve
23: @*/
24: int TaoDAAppSolve(TAO_APPLICATION daapplication, TAO_SOLVER tao){
25: int info;
26: PetscInt i,j,iter;
27: PetscInt mx,my,mz;
28: double fval,residual;
29: TaoTerminateReason reason;
30: DA_APPLICATION daapp;
35: info=TaoAppDAApp(daapplication,&daapp); CHKERRQ(info);
37: if (Tao_DA_Solve==0){
38: info=PetscLogEventRegister("TaoSolveDAApp",DAAPP_COOKIE,&Tao_DA_Solve); CHKERRQ(info);
39: }
41: info = PetscLogEventBegin(Tao_DA_Solve,tao,daapp,0,0);CHKERRQ(info);
42: for (i=0;i<daapp->nda; i++){
43:
44:
45: info=TaoAppResetCounters(daapplication);CHKERRQ(info);
46: info = DAGetInfo(daapp->grid[i].da,PETSC_NULL,&mx,&my,&mz,PETSC_NULL,
47: PETSC_NULL,PETSC_NULL,PETSC_NULL,PETSC_NULL,
48: PETSC_NULL,PETSC_NULL);CHKERRQ(info);
50: info = PetscInfo2(daapp->grid[i].da,"Level %d of %d in TAO_DA_APPLICATION object.\n",i,daapp->nda); CHKERRQ(info);
51: info = PetscInfo3(daapp->grid[i].da,"Global dimensions of DA: mx=%d, my=%d, mz=%d .\n",mx,my,mz); CHKERRQ(info);
53: daapp->currentlevel=i;
55: info=TaoAppSetColoring(daapplication, daapp->grid[i].coloring);CHKERRQ(info);
57: if (i>0){
58: info = TaoSetDown(tao); CHKERRQ(info);
59: info = MatMult(daapp->grid[i].Interpolate,daapp->grid[i-1].X,daapp->grid[i].X); CHKERRQ(info);
60: }
62: info = TaoAppSetInitialSolutionVec(daapplication,daapp->grid[i].X);CHKERRQ(info);
64: if (daapp->grid[i].H){
65: if (!daapplication->ksp) {
66: MPI_Comm comm;
67: info = PetscObjectGetComm((PetscObject)daapp->grid[i].H,&comm); CHKERRQ(info);
68: info = KSPCreate(comm,&daapp->grid[i].ksp); CHKERRQ(info);
69: info = KSPSetFromOptions(daapp->grid[i].ksp); CHKERRQ(info);
70: daapplication->ksp = daapp->grid[i].ksp;
71: /*
72: if (tao->ksp) {
73: tao->ksp->SetKSP(daapp->grid[i].ksp);
74: }
75: info = TaoWrapKSP(daapp->grid[i].ksp,&tao->ksp); CHKERRQ(info);
76: */
77: }
78: if (daapp->IsComplementarity==PETSC_FALSE){
79: info = TaoAppSetHessianMat(daapplication,daapp->grid[i].H,daapp->grid[i].H);CHKERRQ(info);
80: } else {
81: info = TaoAppSetJacobianMat(daapplication,daapp->grid[i].H,daapp->grid[i].H);CHKERRQ(info);
82: }
83: }
85: info = TaoSetupApplicationSolver(daapplication, tao); CHKERRQ(info);
86: for (j=0;j<daapp->nbeforemonitors;j++){
87: info = PetscInfo(daapp->grid[i].da,"TAO: Call before user monitor for DA_APPLICATION object.\n"); CHKERRQ(info);
88: info = (*daapp->beforemonitor[j])(daapplication,daapp->grid[i].da,i,
89: daapp->beforemonitorctx[j]); CHKERRQ(info);
90: }
92: // info = TaoSetUp(tao); CHKERRQ(info);
94: info = PetscInfo1(daapp->grid[i].da,"TAO: Begin solving level %d of DA_APPLICATION object.\n",i); CHKERRQ(info);
95: info = TaoSolveApplication(daapplication,tao);CHKERRQ(info);
97: info = TaoGetSolutionStatus(tao,&iter,&fval,&residual,0,0,&reason); CHKERRQ(info);
98: if (reason <= 0 ){
99: info = PetscInfo1(daapp->grid[i].da,"FAILURE TO FIND SOLUTION at level %d of DA_APPLICATION.\n",i+1); CHKERRQ(info);
100: info = PetscInfo1(daapp->grid[i].da," TAO Reason for termination: %d.\n",reason); CHKERRQ(info);
101: } else {
102: info = PetscInfo1(daapp->grid[i].da,"Found solution to DA_APPLICATION at level: %d.\n",i+1);CHKERRQ(info);
103: info = PetscInfo3(daapp->grid[i].da,"Iterations: %d, Objective Value: %10.8e, Residual: %4.2e.\n",
104: iter,fval,residual);CHKERRQ(info);
105: }
107: for (j=0;j<daapp->naftermonitors;j++){
108: info = PetscInfo(daapp->grid[i].da,"TAO: Call after user monitor for DA_APPLICATION object.\n"); CHKERRQ(info);
109: info = (*daapp->aftermonitor[j])(daapplication,daapp->grid[i].da,i,
110: daapp->aftermonitorctx[j]); CHKERRQ(info);
111: }
113: if (daapplication->ksp) {
114: info = KSPDestroy(daapplication->ksp); CHKERRQ(info);
115: daapplication->ksp=0;
116: }
117: if (i < daapp->nda-1){
118: info = TaoSetDown(tao); CHKERRQ(info);
119: }
120:
121: }
123: info = PetscLogEventEnd(Tao_DA_Solve,tao,daapp,0,0);CHKERRQ(info);
124: return(0);
125: }
127: