1! 2! Description: This example solves a nonlinear system on 1 processor with SNES. 3! We solve the Bratu (SFI - solid fuel ignition) problem in a 2D rectangular 4! domain. The command line options include: 5! -par <parameter>, where <parameter> indicates the nonlinearity of the problem 6! problem SFI: <parameter> = Bratu parameter (0 <= par <= 6.81) 7! -mx <xg>, where <xg> = number of grid points in the x-direction 8! -my <yg>, where <yg> = number of grid points in the y-direction 9! 10 11! 12! -------------------------------------------------------------------------- 13! 14! Solid Fuel Ignition (SFI) problem. This problem is modeled by 15! the partial differential equation 16! 17! -Laplacian u - lambda*exp(u) = 0, 0 < x,y < 1, 18! 19! with boundary conditions 20! 21! u = 0 for x = 0, x = 1, y = 0, y = 1. 22! 23! A finite difference approximation with the usual 5-point stencil 24! is used to discretize the boundary value problem to obtain a nonlinear 25! system of equations. 26! 27! The parallel version of this code is snes/tutorials/ex5f.F 28! 29! -------------------------------------------------------------------------- 30 subroutine postcheck(snes,x,y,w,changed_y,changed_w,ctx,ierr) 31#include <petsc/finclude/petscsnes.h> 32 use petscsnes 33 implicit none 34 SNES snes 35 PetscReal norm 36 Vec tmp,x,y,w 37 PetscBool changed_w,changed_y 38 PetscErrorCode ierr 39 PetscInt ctx 40 PetscScalar mone 41 MPI_Comm comm 42 43 character(len=PETSC_MAX_PATH_LEN) :: outputString 44 45 PetscCallA(PetscObjectGetComm(snes,comm,ierr)) 46 PetscCallA(VecDuplicate(x,tmp,ierr)) 47 mone = -1.0 48 PetscCallA(VecWAXPY(tmp,mone,x,w,ierr)) 49 PetscCallA(VecNorm(tmp,NORM_2,norm,ierr)) 50 PetscCallA(VecDestroy(tmp,ierr)) 51 write(outputString,*) norm 52 PetscCallA(PetscPrintf(comm,'Norm of search step '//trim(outputString)//'\n',ierr)) 53 end 54 55 program main 56#include <petsc/finclude/petscdraw.h> 57 use petscsnes 58 implicit none 59 interface SNESSetJacobian 60 subroutine SNESSetJacobian1(a,b,c,d,e,z) 61 use petscsnes 62 SNES a 63 Mat b 64 Mat c 65 external d 66 MatFDColoring e 67 PetscErrorCode z 68 end subroutine 69 subroutine SNESSetJacobian2(a,b,c,d,e,z) 70 use petscsnes 71 SNES a 72 Mat b 73 Mat c 74 external d 75 integer e 76 PetscErrorCode z 77 end subroutine 78 end interface 79! 80! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 81! Variable declarations 82! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 83! 84! Variables: 85! snes - nonlinear solver 86! x, r - solution, residual vectors 87! J - Jacobian matrix 88! its - iterations for convergence 89! matrix_free - flag - 1 indicates matrix-free version 90! lambda - nonlinearity parameter 91! draw - drawing context 92! 93 SNES snes 94 MatColoring mc 95 Vec x,r 96 PetscDraw draw 97 Mat J 98 PetscBool matrix_free,flg,fd_coloring 99 PetscErrorCode ierr 100 PetscInt its,N, mx,my,i5 101 PetscMPIInt size,rank 102 PetscReal lambda_max,lambda_min,lambda 103 MatFDColoring fdcoloring 104 ISColoring iscoloring 105 PetscBool pc 106 external postcheck 107 108 character(len=PETSC_MAX_PATH_LEN) :: outputString 109 110 PetscScalar,pointer :: lx_v(:) 111 112! Store parameters in common block 113 114 common /params/ lambda,mx,my,fd_coloring 115 116! Note: Any user-defined Fortran routines (such as FormJacobian) 117! MUST be declared as external. 118 119 external FormFunction,FormInitialGuess,FormJacobian 120 121! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 122! Initialize program 123! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 124 125 PetscCallA(PetscInitialize(ierr)) 126 PetscCallMPIA(MPI_Comm_size(PETSC_COMM_WORLD,size,ierr)) 127 PetscCallMPIA(MPI_Comm_rank(PETSC_COMM_WORLD,rank,ierr)) 128 129 PetscCheckA(size .eq. 1,PETSC_COMM_SELF,PETSC_ERR_WRONG_MPI_SIZE,'This is a uniprocessor example only') 130 131! Initialize problem parameters 132 i5 = 5 133 lambda_max = 6.81 134 lambda_min = 0.0 135 lambda = 6.0 136 mx = 4 137 my = 4 138 PetscCallA(PetscOptionsGetInt(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-mx',mx,flg,ierr)) 139 PetscCallA(PetscOptionsGetInt(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-my',my,flg,ierr)) 140 PetscCallA(PetscOptionsGetReal(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-par',lambda,flg,ierr)) 141 PetscCheckA(lambda .lt. lambda_max .and. lambda .gt. lambda_min,PETSC_COMM_SELF,PETSC_ERR_USER,'Lambda out of range ') 142 N = mx*my 143 pc = PETSC_FALSE 144 PetscCallA(PetscOptionsGetBool(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-pc',pc,PETSC_NULL_BOOL,ierr)) 145 146! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 147! Create nonlinear solver context 148! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 149 150 PetscCallA(SNESCreate(PETSC_COMM_WORLD,snes,ierr)) 151 152 if (pc .eqv. PETSC_TRUE) then 153 PetscCallA(SNESSetType(snes,SNESNEWTONTR,ierr)) 154 PetscCallA(SNESNewtonTRSetPostCheck(snes, postcheck,snes,ierr)) 155 endif 156 157! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 158! Create vector data structures; set function evaluation routine 159! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 160 161 PetscCallA(VecCreate(PETSC_COMM_WORLD,x,ierr)) 162 PetscCallA(VecSetSizes(x,PETSC_DECIDE,N,ierr)) 163 PetscCallA(VecSetFromOptions(x,ierr)) 164 PetscCallA(VecDuplicate(x,r,ierr)) 165 166! Set function evaluation routine and vector. Whenever the nonlinear 167! solver needs to evaluate the nonlinear function, it will call this 168! routine. 169! - Note that the final routine argument is the user-defined 170! context that provides application-specific data for the 171! function evaluation routine. 172 173 PetscCallA(SNESSetFunction(snes,r,FormFunction,fdcoloring,ierr)) 174 175! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 176! Create matrix data structure; set Jacobian evaluation routine 177! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 178 179! Create matrix. Here we only approximately preallocate storage space 180! for the Jacobian. See the users manual for a discussion of better 181! techniques for preallocating matrix memory. 182 183 PetscCallA(PetscOptionsHasName(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-snes_mf',matrix_free,ierr)) 184 if (.not. matrix_free) then 185 PetscCallA(MatCreateSeqAIJ(PETSC_COMM_WORLD,N,N,i5,PETSC_NULL_INTEGER,J,ierr)) 186 endif 187 188! 189! This option will cause the Jacobian to be computed via finite differences 190! efficiently using a coloring of the columns of the matrix. 191! 192 fd_coloring = .false. 193 PetscCallA(PetscOptionsHasName(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-snes_fd_coloring',fd_coloring,ierr)) 194 if (fd_coloring) then 195 196! 197! This initializes the nonzero structure of the Jacobian. This is artificial 198! because clearly if we had a routine to compute the Jacobian we won't need 199! to use finite differences. 200! 201 PetscCallA(FormJacobian(snes,x,J,J,0,ierr)) 202! 203! Color the matrix, i.e. determine groups of columns that share no common 204! rows. These columns in the Jacobian can all be computed simultaneously. 205! 206 PetscCallA(MatColoringCreate(J,mc,ierr)) 207 PetscCallA(MatColoringSetType(mc,MATCOLORINGNATURAL,ierr)) 208 PetscCallA(MatColoringSetFromOptions(mc,ierr)) 209 PetscCallA(MatColoringApply(mc,iscoloring,ierr)) 210 PetscCallA(MatColoringDestroy(mc,ierr)) 211! 212! Create the data structure that SNESComputeJacobianDefaultColor() uses 213! to compute the actual Jacobians via finite differences. 214! 215 PetscCallA(MatFDColoringCreate(J,iscoloring,fdcoloring,ierr)) 216 PetscCallA(MatFDColoringSetFunction(fdcoloring,FormFunction,fdcoloring,ierr)) 217 PetscCallA(MatFDColoringSetFromOptions(fdcoloring,ierr)) 218 PetscCallA(MatFDColoringSetUp(J,iscoloring,fdcoloring,ierr)) 219! 220! Tell SNES to use the routine SNESComputeJacobianDefaultColor() 221! to compute Jacobians. 222! 223 PetscCallA(SNESSetJacobian(snes,J,J,SNESComputeJacobianDefaultColor,fdcoloring,ierr)) 224 PetscCallA(ISColoringDestroy(iscoloring,ierr)) 225 226 else if (.not. matrix_free) then 227 228! Set Jacobian matrix data structure and default Jacobian evaluation 229! routine. Whenever the nonlinear solver needs to compute the 230! Jacobian matrix, it will call this routine. 231! - Note that the final routine argument is the user-defined 232! context that provides application-specific data for the 233! Jacobian evaluation routine. 234! - The user can override with: 235! -snes_fd : default finite differencing approximation of Jacobian 236! -snes_mf : matrix-free Newton-Krylov method with no preconditioning 237! (unless user explicitly sets preconditioner) 238! -snes_mf_operator : form preconditioning matrix as set by the user, 239! but use matrix-free approx for Jacobian-vector 240! products within Newton-Krylov method 241! 242 PetscCallA(SNESSetJacobian(snes,J,J,FormJacobian,0,ierr)) 243 endif 244 245! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 246! Customize nonlinear solver; set runtime options 247! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 248 249! Set runtime options (e.g., -snes_monitor -snes_rtol <rtol> -ksp_type <type>) 250 251 PetscCallA(SNESSetFromOptions(snes,ierr)) 252 253! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 254! Evaluate initial guess; then solve nonlinear system. 255! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 256 257! Note: The user should initialize the vector, x, with the initial guess 258! for the nonlinear solver prior to calling SNESSolve(). In particular, 259! to employ an initial guess of zero, the user should explicitly set 260! this vector to zero by calling VecSet(). 261 262 PetscCallA(FormInitialGuess(x,ierr)) 263 PetscCallA(SNESSolve(snes,PETSC_NULL_VEC,x,ierr)) 264 PetscCallA(SNESGetIterationNumber(snes,its,ierr)) 265 write(outputString,*) its 266 PetscCallA(PetscPrintf(PETSC_COMM_WORLD,'Number of SNES iterations = '//trim(outputString)//'\n',ierr)) 267 268! PetscDraw contour plot of solution 269 270 PetscCallA(PetscDrawCreate(PETSC_COMM_WORLD,PETSC_NULL_CHARACTER,'Solution',300,0,300,300,draw,ierr)) 271 PetscCallA(PetscDrawSetFromOptions(draw,ierr)) 272 273 PetscCallA(VecGetArrayReadF90(x,lx_v,ierr)) 274 PetscCallA(PetscDrawTensorContour(draw,mx,my,PETSC_NULL_REAL,PETSC_NULL_REAL,lx_v,ierr)) 275 PetscCallA(VecRestoreArrayReadF90(x,lx_v,ierr)) 276 277! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 278! Free work space. All PETSc objects should be destroyed when they 279! are no longer needed. 280! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 281 282 if (.not. matrix_free) PetscCallA(MatDestroy(J,ierr)) 283 if (fd_coloring) PetscCallA(MatFDColoringDestroy(fdcoloring,ierr)) 284 285 PetscCallA(VecDestroy(x,ierr)) 286 PetscCallA(VecDestroy(r,ierr)) 287 PetscCallA(SNESDestroy(snes,ierr)) 288 PetscCallA(PetscDrawDestroy(draw,ierr)) 289 PetscCallA(PetscFinalize(ierr)) 290 end 291 292! --------------------------------------------------------------------- 293! 294! FormInitialGuess - Forms initial approximation. 295! 296! Input Parameter: 297! X - vector 298! 299! Output Parameters: 300! X - vector 301! ierr - error code 302! 303! Notes: 304! This routine serves as a wrapper for the lower-level routine 305! "ApplicationInitialGuess", where the actual computations are 306! done using the standard Fortran style of treating the local 307! vector data as a multidimensional array over the local mesh. 308! This routine merely accesses the local vector data via 309! VecGetArrayF90() and VecRestoreArrayF90(). 310! 311 subroutine FormInitialGuess(X,ierr) 312 use petscsnes 313 implicit none 314 315! Input/output variables: 316 Vec X 317 PetscErrorCode ierr 318 319! Declarations for use with local arrays: 320 PetscScalar,pointer :: lx_v(:) 321 322 ierr = 0 323 324! Get a pointer to vector data. 325! - VecGetArrayF90() returns a pointer to the data array. 326! - You MUST call VecRestoreArrayF90() when you no longer need access to 327! the array. 328 329 PetscCallA(VecGetArrayF90(X,lx_v,ierr)) 330 331! Compute initial guess 332 333 PetscCallA(ApplicationInitialGuess(lx_v,ierr)) 334 335! Restore vector 336 337 PetscCallA(VecRestoreArrayF90(X,lx_v,ierr)) 338 339 end 340 341! ApplicationInitialGuess - Computes initial approximation, called by 342! the higher level routine FormInitialGuess(). 343! 344! Input Parameter: 345! x - local vector data 346! 347! Output Parameters: 348! f - local vector data, f(x) 349! ierr - error code 350! 351! Notes: 352! This routine uses standard Fortran-style computations over a 2-dim array. 353! 354 subroutine ApplicationInitialGuess(x,ierr) 355 use petscksp 356 implicit none 357 358! Common blocks: 359 PetscReal lambda 360 PetscInt mx,my 361 PetscBool fd_coloring 362 common /params/ lambda,mx,my,fd_coloring 363 364! Input/output variables: 365 PetscScalar x(mx,my) 366 PetscErrorCode ierr 367 368! Local variables: 369 PetscInt i,j 370 PetscReal temp1,temp,hx,hy,one 371 372! Set parameters 373 374 ierr = 0 375 one = 1.0 376 hx = one/(mx-1) 377 hy = one/(my-1) 378 temp1 = lambda/(lambda + one) 379 380 do 20 j=1,my 381 temp = min(j-1,my-j)*hy 382 do 10 i=1,mx 383 if (i .eq. 1 .or. j .eq. 1 .or. i .eq. mx .or. j .eq. my) then 384 x(i,j) = 0.0 385 else 386 x(i,j) = temp1 * sqrt(min(min(i-1,mx-i)*hx,temp)) 387 endif 388 10 continue 389 20 continue 390 391 end 392 393! --------------------------------------------------------------------- 394! 395! FormFunction - Evaluates nonlinear function, F(x). 396! 397! Input Parameters: 398! snes - the SNES context 399! X - input vector 400! dummy - optional user-defined context, as set by SNESSetFunction() 401! (not used here) 402! 403! Output Parameter: 404! F - vector with newly computed function 405! 406! Notes: 407! This routine serves as a wrapper for the lower-level routine 408! "ApplicationFunction", where the actual computations are 409! done using the standard Fortran style of treating the local 410! vector data as a multidimensional array over the local mesh. 411! This routine merely accesses the local vector data via 412! VecGetArrayF90() and VecRestoreArrayF90(). 413! 414 subroutine FormFunction(snes,X,F,fdcoloring,ierr) 415 use petscsnes 416 implicit none 417 418! Input/output variables: 419 SNES snes 420 Vec X,F 421 PetscErrorCode ierr 422 MatFDColoring fdcoloring 423 424! Common blocks: 425 PetscReal lambda 426 PetscInt mx,my 427 PetscBool fd_coloring 428 common /params/ lambda,mx,my,fd_coloring 429 430! Declarations for use with local arrays: 431 PetscScalar,pointer :: lx_v(:), lf_v(:) 432 PetscInt, pointer :: indices(:) 433 434! Get pointers to vector data. 435! - VecGetArrayF90() returns a pointer to the data array. 436! - You MUST call VecRestoreArrayF90() when you no longer need access to 437! the array. 438 439 PetscCallA(VecGetArrayReadF90(X,lx_v,ierr)) 440 PetscCallA(VecGetArrayF90(F,lf_v,ierr)) 441 442! Compute function 443 444 PetscCallA(ApplicationFunction(lx_v,lf_v,ierr)) 445 446! Restore vectors 447 448 PetscCallA(VecRestoreArrayReadF90(X,lx_v,ierr)) 449 PetscCallA(VecRestoreArrayF90(F,lf_v,ierr)) 450 451 PetscCallA(PetscLogFlops(11.0d0*mx*my,ierr)) 452! 453! fdcoloring is in the common block and used here ONLY to test the 454! calls to MatFDColoringGetPerturbedColumnsF90() and MatFDColoringRestorePerturbedColumnsF90() 455! 456 if (fd_coloring) then 457 PetscCallA(MatFDColoringGetPerturbedColumnsF90(fdcoloring,indices,ierr)) 458 print*,'Indices from GetPerturbedColumnsF90' 459 write(*,1000) indices 460 1000 format(50i4) 461 PetscCallA(MatFDColoringRestorePerturbedColumnsF90(fdcoloring,indices,ierr)) 462 endif 463 end 464 465! --------------------------------------------------------------------- 466! 467! ApplicationFunction - Computes nonlinear function, called by 468! the higher level routine FormFunction(). 469! 470! Input Parameter: 471! x - local vector data 472! 473! Output Parameters: 474! f - local vector data, f(x) 475! ierr - error code 476! 477! Notes: 478! This routine uses standard Fortran-style computations over a 2-dim array. 479! 480 subroutine ApplicationFunction(x,f,ierr) 481 use petscsnes 482 implicit none 483 484! Common blocks: 485 PetscReal lambda 486 PetscInt mx,my 487 PetscBool fd_coloring 488 common /params/ lambda,mx,my,fd_coloring 489 490! Input/output variables: 491 PetscScalar x(mx,my),f(mx,my) 492 PetscErrorCode ierr 493 494! Local variables: 495 PetscScalar two,one,hx,hy 496 PetscScalar hxdhy,hydhx,sc 497 PetscScalar u,uxx,uyy 498 PetscInt i,j 499 500 ierr = 0 501 one = 1.0 502 two = 2.0 503 hx = one/(mx-1) 504 hy = one/(my-1) 505 sc = hx*hy*lambda 506 hxdhy = hx/hy 507 hydhx = hy/hx 508 509! Compute function 510 511 do 20 j=1,my 512 do 10 i=1,mx 513 if (i .eq. 1 .or. j .eq. 1 .or. i .eq. mx .or. j .eq. my) then 514 f(i,j) = x(i,j) 515 else 516 u = x(i,j) 517 uxx = hydhx * (two*u - x(i-1,j) - x(i+1,j)) 518 uyy = hxdhy * (two*u - x(i,j-1) - x(i,j+1)) 519 f(i,j) = uxx + uyy - sc*exp(u) 520 endif 521 10 continue 522 20 continue 523 524 end 525 526! --------------------------------------------------------------------- 527! 528! FormJacobian - Evaluates Jacobian matrix. 529! 530! Input Parameters: 531! snes - the SNES context 532! x - input vector 533! dummy - optional user-defined context, as set by SNESSetJacobian() 534! (not used here) 535! 536! Output Parameters: 537! jac - Jacobian matrix 538! jac_prec - optionally different preconditioning matrix (not used here) 539! flag - flag indicating matrix structure 540! 541! Notes: 542! This routine serves as a wrapper for the lower-level routine 543! "ApplicationJacobian", where the actual computations are 544! done using the standard Fortran style of treating the local 545! vector data as a multidimensional array over the local mesh. 546! This routine merely accesses the local vector data via 547! VecGetArrayF90() and VecRestoreArrayF90(). 548! 549 subroutine FormJacobian(snes,X,jac,jac_prec,dummy,ierr) 550 use petscsnes 551 implicit none 552 553! Input/output variables: 554 SNES snes 555 Vec X 556 Mat jac,jac_prec 557 PetscErrorCode ierr 558 integer dummy 559 560! Common blocks: 561 PetscReal lambda 562 PetscInt mx,my 563 PetscBool fd_coloring 564 common /params/ lambda,mx,my,fd_coloring 565 566! Declarations for use with local array: 567 PetscScalar,pointer :: lx_v(:) 568 569! Get a pointer to vector data 570 571 PetscCallA(VecGetArrayReadF90(X,lx_v,ierr)) 572 573! Compute Jacobian entries 574 575 PetscCallA(ApplicationJacobian(lx_v,jac,jac_prec,ierr)) 576 577! Restore vector 578 579 PetscCallA(VecRestoreArrayReadF90(X,lx_v,ierr)) 580 581! Assemble matrix 582 583 PetscCallA(MatAssemblyBegin(jac_prec,MAT_FINAL_ASSEMBLY,ierr)) 584 PetscCallA(MatAssemblyEnd(jac_prec,MAT_FINAL_ASSEMBLY,ierr)) 585 586 end 587 588! --------------------------------------------------------------------- 589! 590! ApplicationJacobian - Computes Jacobian matrix, called by 591! the higher level routine FormJacobian(). 592! 593! Input Parameters: 594! x - local vector data 595! 596! Output Parameters: 597! jac - Jacobian matrix 598! jac_prec - optionally different preconditioning matrix (not used here) 599! ierr - error code 600! 601! Notes: 602! This routine uses standard Fortran-style computations over a 2-dim array. 603! 604 subroutine ApplicationJacobian(x,jac,jac_prec,ierr) 605 use petscsnes 606 implicit none 607 608! Common blocks: 609 PetscReal lambda 610 PetscInt mx,my 611 PetscBool fd_coloring 612 common /params/ lambda,mx,my,fd_coloring 613 614! Input/output variables: 615 PetscScalar x(mx,my) 616 Mat jac,jac_prec 617 PetscErrorCode ierr 618 619! Local variables: 620 PetscInt i,j,row(1),col(5),i1,i5 621 PetscScalar two,one, hx,hy 622 PetscScalar hxdhy,hydhx,sc,v(5) 623 624! Set parameters 625 626 i1 = 1 627 i5 = 5 628 one = 1.0 629 two = 2.0 630 hx = one/(mx-1) 631 hy = one/(my-1) 632 sc = hx*hy 633 hxdhy = hx/hy 634 hydhx = hy/hx 635 636! Compute entries of the Jacobian matrix 637! - Here, we set all entries for a particular row at once. 638! - Note that MatSetValues() uses 0-based row and column numbers 639! in Fortran as well as in C. 640 641 do 20 j=1,my 642 row(1) = (j-1)*mx - 1 643 do 10 i=1,mx 644 row(1) = row(1) + 1 645! boundary points 646 if (i .eq. 1 .or. j .eq. 1 .or. i .eq. mx .or. j .eq. my) then 647 PetscCallA(MatSetValues(jac_prec,i1,row,i1,row,one,INSERT_VALUES,ierr)) 648! interior grid points 649 else 650 v(1) = -hxdhy 651 v(2) = -hydhx 652 v(3) = two*(hydhx + hxdhy) - sc*lambda*exp(x(i,j)) 653 v(4) = -hydhx 654 v(5) = -hxdhy 655 col(1) = row(1) - mx 656 col(2) = row(1) - 1 657 col(3) = row(1) 658 col(4) = row(1) + 1 659 col(5) = row(1) + mx 660 PetscCallA(MatSetValues(jac_prec,i1,row,i5,col,v,INSERT_VALUES,ierr)) 661 endif 662 10 continue 663 20 continue 664 665 end 666 667! 668!/*TEST 669! 670! build: 671! requires: !single 672! 673! test: 674! args: -snes_monitor_short -nox -snes_type newtontr -ksp_gmres_cgs_refinement_type refine_always 675! 676! test: 677! suffix: 2 678! args: -snes_monitor_short -nox -snes_fd -ksp_gmres_cgs_refinement_type refine_always 679! 680! test: 681! suffix: 3 682! args: -snes_monitor_short -nox -snes_fd_coloring -mat_coloring_type sl -ksp_gmres_cgs_refinement_type refine_always 683! filter: sort -b 684! filter_output: sort -b 685! 686! test: 687! suffix: 4 688! args: -pc -par 6.807 -nox 689! 690!TEST*/ 691