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