Created
March 10, 2015 14:40
-
-
Save jhale/f97256a9aab4a54430ba to your computer and use it in GitHub Desktop.
Example of how to pass command scoped options to PETSc and SLEPc using DOLFIN
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| import sys | |
| # Here I make two namespaces, ns_ and A_. | |
| # You can make as many as you want, probably | |
| # one for every KSP/SVD/EPS in your code. | |
| args = [sys.argv[0]] + """ | |
| --petsc.ns_svd_type lapack | |
| --petsc.A_mat_type aij | |
| """.split() | |
| # Don't init petsc4py and slepc4py yourself, let DOLFIN handle it. | |
| # You will run into problems with MPI communicators if you do it yourself. | |
| from dolfin import * | |
| # This line is also crucial! | |
| parameters.parse(args) | |
| # DOLFIN has already run init on these for you. | |
| from petsc4py import PETSc | |
| from slepc4py import SLEPc | |
| n = 60 | |
| mu = 1e-06 | |
| print("Lauching singular value decomposition, (%d x %d) mu=%g\n" % (n+1, n, mu) ) | |
| A = PETSc.Mat(); A.create() | |
| A.setSizes([n+1, n]) | |
| # Here I set the options prefix, remember the trailing _! | |
| A.setOptionsPrefix("A_") | |
| # Here the options are parsed from the A_ prefix. | |
| A.setFromOptions() | |
| A.setUp() | |
| rstart, rend = A.getOwnershipRange() | |
| for i in xrange(rstart, rend): | |
| if i==0: | |
| for j in range(n): | |
| A[0,j] = 1.0 | |
| else: | |
| A[i,i-1] = mu | |
| A.assemble() | |
| S = SLEPc.SVD() | |
| S.create() | |
| S.setOperator(A) | |
| # and here I set the options prefix. | |
| S.setOptionsPrefix("ns_") | |
| S.setFromOptions() | |
| S.solve() | |
| Print = PETSc.Sys.Print | |
| Print( "******************************" ) | |
| Print( "*** SLEPc Solution Results ***" ) | |
| Print( "******************************\n" ) | |
| svd_type = S.getType() | |
| Print( "Solution method: %s" % svd_type ) | |
| its = S.getIterationNumber() | |
| Print( "Number of iterations of the method: %d" % its ) | |
| nsv, ncv, mpd = S.getDimensions() | |
| Print( "Number of requested singular values: %d" % nsv ) | |
| tol, maxit = S.getTolerances() | |
| Print( "Stopping condition: tol=%.4g, maxit=%d" % (tol, maxit) ) | |
| nconv = S.getConverged() | |
| Print( "Number of converged approximate singular triplets %d" % nconv ) | |
| if nconv > 0: | |
| v, u = A.getVecs() | |
| Print() | |
| Print(" sigma residual norm ") | |
| Print("------------- ---------------") | |
| for i in range(nconv): | |
| sigma = S.getSingularTriplet(i, u, v) | |
| error = S.computeRelativeError(i) | |
| Print( " %6f %12g" % (sigma, error) ) | |
| Print() |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment