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 "daapp_impl.h"
7: static int Tao_DA_Solve=0;
9: int TaoAppDAApp(TAO_APPLICATION, DA_APPLICATION *);
13: /*@
14: TaoDAAppSolve - Solve the PETSC DA application.
16: Input Parameters:
17: . daapplication - the TAO DAApplication structure
19: Level: beginner
21: .keywords: Application, DA, Solve
22: @*/
23: int TaoDAAppSolve(TAO_APPLICATION daapplication, TAO_SOLVER tao){
24: int i,j,mx,my,mz,iter,info;
25: double fval,residual;
26: TaoTerminateReason reason;
27: DA_APPLICATION daapp;
32: info=TaoAppDAApp(daapplication,&daapp); CHKERRQ(info);
34: if (Tao_DA_Solve==0){
35: info=PetscLogEventRegister(&Tao_DA_Solve,"TaoSolveDAApp",DAAPP_COOKIE); CHKERRQ(info);
36: }
38: info = PetscLogEventBegin(Tao_DA_Solve,tao,daapp,0,0);CHKERRQ(info);
39: for (i=0;i<daapp->nda; i++){
40:
41: info=TaoAppResetCounters(daapplication);CHKERRQ(info);
42: info = DAGetInfo(daapp->grid[i].da,PETSC_NULL,&mx,&my,&mz,PETSC_NULL,
43: PETSC_NULL,PETSC_NULL,PETSC_NULL,PETSC_NULL,
44: PETSC_NULL,PETSC_NULL);CHKERRQ(info);
46: info = PetscLogInfo((daapp->grid[i].da,"Level %d of %d in TAO_DA_APPLICATION object.\n",i,daapp->nda)); CHKERRQ(info);
47: info = PetscLogInfo((daapp->grid[i].da,"Global dimensions of DA: mx=%d, my=%d, mz=%d .\n",mx,my,mz)); CHKERRQ(info);
49: daapp->currentlevel=i;
51: info=TaoAppSetColoring(daapplication, daapp->grid[i].coloring);CHKERRQ(info);
53: if (i>0){
54: info = TaoSetDown(tao); CHKERRQ(info);
55: info = MatMult(daapp->grid[i].Interpolate,daapp->grid[i-1].X,daapp->grid[i].X); CHKERRQ(info);
56: }
58: info = TaoAppSetInitialSolutionVec(daapplication,daapp->grid[i].X);CHKERRQ(info);
60: if (daapp->grid[i].H){
61: if (daapp->IsComplementarity==PETSC_FALSE){
62: info = TaoAppSetHessianMat(daapplication,daapp->grid[i].H,daapp->grid[i].H);CHKERRQ(info);
63: } else {
64: info = TaoAppSetJacobianMat(daapplication,daapp->grid[i].H,daapp->grid[i].H);CHKERRQ(info);
65: }
66: }
68: info = TaoSetupApplicationSolver(daapplication, tao); CHKERRQ(info);
69: for (j=0;j<daapp->nbeforemonitors;j++){
70: info = PetscLogInfo((daapp->grid[i].da,"TAO: Call before user monitor for DA_APPLICATION object.\n")); CHKERRQ(info);
71: info = (*daapp->beforemonitor[j])(daapplication,daapp->grid[i].da,i,
72: daapp->beforemonitorctx[j]); CHKERRQ(info);
73: }
75: // info = TaoSetUp(tao); CHKERRQ(info);
77: info = PetscLogInfo((daapp->grid[i].da,"TAO: Begin solving level %d of DA_APPLICATION object.\n",i)); CHKERRQ(info);
78: info = TaoSolveApplication(daapplication,tao);CHKERRQ(info);
80: info = TaoGetSolutionStatus(tao,&iter,&fval,&residual,0,0,&reason); CHKERRQ(info);
81: if (reason <= 0 ){
82: info = PetscLogInfo((daapp->grid[i].da,"FAILURE TO FIND SOLUTION at level %d of DA_APPLICATION.\n",i+1)); CHKERRQ(info);
83: info = PetscLogInfo((daapp->grid[i].da," TAO Reason for termination: %d.\n",reason)); CHKERRQ(info);
84: } else {
85: info = PetscLogInfo((daapp->grid[i].da,"Found solution to DA_APPLICATION at level: %d.\n",i+1));CHKERRQ(info);
86: info = PetscLogInfo((daapp->grid[i].da,"Iterations: %d, Objective Value: %10.8e, Residual: %4.2e.\n",
87: iter,fval,residual));CHKERRQ(info);
88: }
90: for (j=0;j<daapp->naftermonitors;j++){
91: info = PetscLogInfo((daapp->grid[i].da,"TAO: Call after user monitor for DA_APPLICATION object.\n")); CHKERRQ(info);
92: info = (*daapp->aftermonitor[j])(daapplication,daapp->grid[i].da,i,
93: daapp->aftermonitorctx[j]); CHKERRQ(info);
94: }
96: if (i < daapp->nda-1){
97: info = TaoSetDown(tao); CHKERRQ(info);
98: }
99:
100: }
102: info = PetscLogEventEnd(Tao_DA_Solve,tao,daapp,0,0);CHKERRQ(info);
103: return(0);
104: }
106: