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