LCOV - code coverage report
Current view: top level - svd/interface - svdbasic.c (source / functions) Hit Total Coverage
Test: SLEPc Lines: 187 198 94.4 %
Date: 2021-08-02 00:32:28 Functions: 10 11 90.9 %
Legend: Lines: hit not hit

          Line data    Source code
       1             : /*
       2             :    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
       3             :    SLEPc - Scalable Library for Eigenvalue Problem Computations
       4             :    Copyright (c) 2002-2021, Universitat Politecnica de Valencia, Spain
       5             : 
       6             :    This file is part of SLEPc.
       7             :    SLEPc is distributed under a 2-clause BSD license (see LICENSE).
       8             :    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
       9             : */
      10             : /*
      11             :    Basic SVD routines
      12             : */
      13             : 
      14             : #include <slepc/private/svdimpl.h>      /*I "slepcsvd.h" I*/
      15             : 
      16             : /* Logging support */
      17             : PetscClassId      SVD_CLASSID = 0;
      18             : PetscLogEvent     SVD_SetUp = 0,SVD_Solve = 0;
      19             : 
      20             : /* List of registered SVD routines */
      21             : PetscFunctionList SVDList = 0;
      22             : PetscBool         SVDRegisterAllCalled = PETSC_FALSE;
      23             : 
      24             : /* List of registered SVD monitors */
      25             : PetscFunctionList SVDMonitorList              = NULL;
      26             : PetscFunctionList SVDMonitorCreateList        = NULL;
      27             : PetscFunctionList SVDMonitorDestroyList       = NULL;
      28             : PetscBool         SVDMonitorRegisterAllCalled = PETSC_FALSE;
      29             : 
      30             : /*@
      31             :    SVDCreate - Creates the default SVD context.
      32             : 
      33             :    Collective
      34             : 
      35             :    Input Parameter:
      36             : .  comm - MPI communicator
      37             : 
      38             :    Output Parameter:
      39             : .  svd - location to put the SVD context
      40             : 
      41             :    Note:
      42             :    The default SVD type is SVDCROSS
      43             : 
      44             :    Level: beginner
      45             : 
      46             : .seealso: SVDSetUp(), SVDSolve(), SVDDestroy(), SVD
      47             : @*/
      48         143 : PetscErrorCode SVDCreate(MPI_Comm comm,SVD *outsvd)
      49             : {
      50         143 :   PetscErrorCode ierr;
      51         143 :   SVD            svd;
      52             : 
      53         143 :   PetscFunctionBegin;
      54         143 :   PetscValidPointer(outsvd,2);
      55         143 :   *outsvd = 0;
      56         143 :   ierr = SVDInitializePackage();CHKERRQ(ierr);
      57         143 :   ierr = SlepcHeaderCreate(svd,SVD_CLASSID,"SVD","Singular Value Decomposition","SVD",comm,SVDDestroy,SVDView);CHKERRQ(ierr);
      58             : 
      59         143 :   svd->OP               = NULL;
      60         143 :   svd->OPb              = NULL;
      61         143 :   svd->max_it           = PETSC_DEFAULT;
      62         143 :   svd->nsv              = 1;
      63         143 :   svd->ncv              = PETSC_DEFAULT;
      64         143 :   svd->mpd              = PETSC_DEFAULT;
      65         143 :   svd->nini             = 0;
      66         143 :   svd->ninil            = 0;
      67         143 :   svd->tol              = PETSC_DEFAULT;
      68         143 :   svd->conv             = SVD_CONV_REL;
      69         143 :   svd->stop             = SVD_STOP_BASIC;
      70         143 :   svd->which            = SVD_LARGEST;
      71         143 :   svd->problem_type     = (SVDProblemType)0;
      72         143 :   svd->impltrans        = PETSC_FALSE;
      73         143 :   svd->trackall         = PETSC_FALSE;
      74             : 
      75         143 :   svd->converged        = SVDConvergedRelative;
      76         143 :   svd->convergeduser    = NULL;
      77         143 :   svd->convergeddestroy = NULL;
      78         143 :   svd->stopping         = SVDStoppingBasic;
      79         143 :   svd->stoppinguser     = NULL;
      80         143 :   svd->stoppingdestroy  = NULL;
      81         143 :   svd->convergedctx     = NULL;
      82         143 :   svd->stoppingctx      = NULL;
      83         143 :   svd->numbermonitors   = 0;
      84             : 
      85         143 :   svd->ds               = NULL;
      86         143 :   svd->U                = NULL;
      87         143 :   svd->V                = NULL;
      88         143 :   svd->A                = NULL;
      89         143 :   svd->B                = NULL;
      90         143 :   svd->AT               = NULL;
      91         143 :   svd->BT               = NULL;
      92         143 :   svd->IS               = NULL;
      93         143 :   svd->ISL              = NULL;
      94         143 :   svd->sigma            = NULL;
      95         143 :   svd->errest           = NULL;
      96         143 :   svd->perm             = NULL;
      97         143 :   svd->nworkl           = 0;
      98         143 :   svd->nworkr           = 0;
      99         143 :   svd->workl            = NULL;
     100         143 :   svd->workr            = NULL;
     101         143 :   svd->data             = NULL;
     102             : 
     103         143 :   svd->state            = SVD_STATE_INITIAL;
     104         143 :   svd->nconv            = 0;
     105         143 :   svd->its              = 0;
     106         143 :   svd->leftbasis        = PETSC_FALSE;
     107         143 :   svd->swapped          = PETSC_FALSE;
     108         143 :   svd->expltrans        = PETSC_FALSE;
     109         143 :   svd->nrma             = 0.0;
     110         143 :   svd->nrmb             = 0.0;
     111         143 :   svd->isgeneralized    = PETSC_FALSE;
     112         143 :   svd->reason           = SVD_CONVERGED_ITERATING;
     113             : 
     114         143 :   ierr = PetscNewLog(svd,&svd->sc);CHKERRQ(ierr);
     115         143 :   *outsvd = svd;
     116         143 :   PetscFunctionReturn(0);
     117             : }
     118             : 
     119             : /*@
     120             :    SVDReset - Resets the SVD context to the initial state (prior to setup)
     121             :    and destroys any allocated Vecs and Mats.
     122             : 
     123             :    Collective on svd
     124             : 
     125             :    Input Parameter:
     126             : .  svd - singular value solver context obtained from SVDCreate()
     127             : 
     128             :    Level: advanced
     129             : 
     130             : .seealso: SVDDestroy()
     131             : @*/
     132         149 : PetscErrorCode SVDReset(SVD svd)
     133             : {
     134         149 :   PetscErrorCode ierr;
     135             : 
     136         149 :   PetscFunctionBegin;
     137         149 :   if (svd) PetscValidHeaderSpecific(svd,SVD_CLASSID,1);
     138         149 :   if (!svd) PetscFunctionReturn(0);
     139         149 :   if (svd->ops->reset) { ierr = (svd->ops->reset)(svd);CHKERRQ(ierr); }
     140         149 :   ierr = MatDestroy(&svd->OP);CHKERRQ(ierr);
     141         149 :   ierr = MatDestroy(&svd->OPb);CHKERRQ(ierr);
     142         149 :   ierr = MatDestroy(&svd->A);CHKERRQ(ierr);
     143         149 :   ierr = MatDestroy(&svd->B);CHKERRQ(ierr);
     144         149 :   ierr = MatDestroy(&svd->AT);CHKERRQ(ierr);
     145         149 :   ierr = MatDestroy(&svd->BT);CHKERRQ(ierr);
     146         149 :   ierr = BVDestroy(&svd->U);CHKERRQ(ierr);
     147         149 :   ierr = BVDestroy(&svd->V);CHKERRQ(ierr);
     148         149 :   ierr = VecDestroyVecs(svd->nworkl,&svd->workl);CHKERRQ(ierr);
     149         149 :   svd->nworkl = 0;
     150         149 :   ierr = VecDestroyVecs(svd->nworkr,&svd->workr);CHKERRQ(ierr);
     151         149 :   svd->nworkr = 0;
     152         149 :   svd->state = SVD_STATE_INITIAL;
     153         149 :   PetscFunctionReturn(0);
     154             : }
     155             : 
     156             : /*@C
     157             :    SVDDestroy - Destroys the SVD context.
     158             : 
     159             :    Collective on svd
     160             : 
     161             :    Input Parameter:
     162             : .  svd - singular value solver context obtained from SVDCreate()
     163             : 
     164             :    Level: beginner
     165             : 
     166             : .seealso: SVDCreate(), SVDSetUp(), SVDSolve()
     167             : @*/
     168         143 : PetscErrorCode SVDDestroy(SVD *svd)
     169             : {
     170         143 :   PetscErrorCode ierr;
     171             : 
     172         143 :   PetscFunctionBegin;
     173         143 :   if (!*svd) PetscFunctionReturn(0);
     174         143 :   PetscValidHeaderSpecific(*svd,SVD_CLASSID,1);
     175         143 :   if (--((PetscObject)(*svd))->refct > 0) { *svd = 0; PetscFunctionReturn(0); }
     176         143 :   ierr = SVDReset(*svd);CHKERRQ(ierr);
     177         143 :   if ((*svd)->ops->destroy) { ierr = (*(*svd)->ops->destroy)(*svd);CHKERRQ(ierr); }
     178         143 :   if ((*svd)->sigma) {
     179         143 :     ierr = PetscFree3((*svd)->sigma,(*svd)->perm,(*svd)->errest);CHKERRQ(ierr);
     180             :   }
     181         143 :   ierr = DSDestroy(&(*svd)->ds);CHKERRQ(ierr);
     182         143 :   ierr = PetscFree((*svd)->sc);CHKERRQ(ierr);
     183             :   /* just in case the initial vectors have not been used */
     184         143 :   ierr = SlepcBasisDestroy_Private(&(*svd)->nini,&(*svd)->IS);CHKERRQ(ierr);
     185         143 :   ierr = SlepcBasisDestroy_Private(&(*svd)->ninil,&(*svd)->ISL);CHKERRQ(ierr);
     186         143 :   ierr = SVDMonitorCancel(*svd);CHKERRQ(ierr);
     187         143 :   ierr = PetscHeaderDestroy(svd);CHKERRQ(ierr);
     188         143 :   PetscFunctionReturn(0);
     189             : }
     190             : 
     191             : /*@C
     192             :    SVDSetType - Selects the particular solver to be used in the SVD object.
     193             : 
     194             :    Logically Collective on svd
     195             : 
     196             :    Input Parameters:
     197             : +  svd      - the singular value solver context
     198             : -  type     - a known method
     199             : 
     200             :    Options Database Key:
     201             : .  -svd_type <method> - Sets the method; use -help for a list
     202             :     of available methods
     203             : 
     204             :    Notes:
     205             :    See "slepc/include/slepcsvd.h" for available methods. The default
     206             :    is SVDCROSS.
     207             : 
     208             :    Normally, it is best to use the SVDSetFromOptions() command and
     209             :    then set the SVD type from the options database rather than by using
     210             :    this routine.  Using the options database provides the user with
     211             :    maximum flexibility in evaluating the different available methods.
     212             :    The SVDSetType() routine is provided for those situations where it
     213             :    is necessary to set the iterative solver independently of the command
     214             :    line or options database.
     215             : 
     216             :    Level: intermediate
     217             : 
     218             : .seealso: SVDType
     219             : @*/
     220         144 : PetscErrorCode SVDSetType(SVD svd,SVDType type)
     221             : {
     222         144 :   PetscErrorCode ierr,(*r)(SVD);
     223         144 :   PetscBool      match;
     224             : 
     225         144 :   PetscFunctionBegin;
     226         144 :   PetscValidHeaderSpecific(svd,SVD_CLASSID,1);
     227         144 :   PetscValidCharPointer(type,2);
     228             : 
     229         144 :   ierr = PetscObjectTypeCompare((PetscObject)svd,type,&match);CHKERRQ(ierr);
     230         144 :   if (match) PetscFunctionReturn(0);
     231             : 
     232         144 :   ierr = PetscFunctionListFind(SVDList,type,&r);CHKERRQ(ierr);
     233         144 :   if (!r) SETERRQ1(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_UNKNOWN_TYPE,"Unknown SVD type given: %s",type);
     234             : 
     235         144 :   if (svd->ops->destroy) { ierr = (*svd->ops->destroy)(svd);CHKERRQ(ierr); }
     236         144 :   ierr = PetscMemzero(svd->ops,sizeof(struct _SVDOps));CHKERRQ(ierr);
     237             : 
     238         144 :   svd->state = SVD_STATE_INITIAL;
     239         144 :   ierr = PetscObjectChangeTypeName((PetscObject)svd,type);CHKERRQ(ierr);
     240         144 :   ierr = (*r)(svd);CHKERRQ(ierr);
     241         144 :   PetscFunctionReturn(0);
     242             : }
     243             : 
     244             : /*@C
     245             :    SVDGetType - Gets the SVD type as a string from the SVD object.
     246             : 
     247             :    Not Collective
     248             : 
     249             :    Input Parameter:
     250             : .  svd - the singular value solver context
     251             : 
     252             :    Output Parameter:
     253             : .  name - name of SVD method
     254             : 
     255             :    Level: intermediate
     256             : 
     257             : .seealso: SVDSetType()
     258             : @*/
     259          35 : PetscErrorCode SVDGetType(SVD svd,SVDType *type)
     260             : {
     261          35 :   PetscFunctionBegin;
     262          35 :   PetscValidHeaderSpecific(svd,SVD_CLASSID,1);
     263          35 :   PetscValidPointer(type,2);
     264          35 :   *type = ((PetscObject)svd)->type_name;
     265          35 :   PetscFunctionReturn(0);
     266             : }
     267             : 
     268             : /*@C
     269             :    SVDRegister - Adds a method to the singular value solver package.
     270             : 
     271             :    Not Collective
     272             : 
     273             :    Input Parameters:
     274             : +  name - name of a new user-defined solver
     275             : -  function - routine to create the solver context
     276             : 
     277             :    Notes:
     278             :    SVDRegister() may be called multiple times to add several user-defined solvers.
     279             : 
     280             :    Sample usage:
     281             : .vb
     282             :     SVDRegister("my_solver",MySolverCreate);
     283             : .ve
     284             : 
     285             :    Then, your solver can be chosen with the procedural interface via
     286             : $     SVDSetType(svd,"my_solver")
     287             :    or at runtime via the option
     288             : $     -svd_type my_solver
     289             : 
     290             :    Level: advanced
     291             : 
     292             : .seealso: SVDRegisterAll()
     293             : @*/
     294        1008 : PetscErrorCode SVDRegister(const char *name,PetscErrorCode (*function)(SVD))
     295             : {
     296        1008 :   PetscErrorCode ierr;
     297             : 
     298        1008 :   PetscFunctionBegin;
     299        1008 :   ierr = SVDInitializePackage();CHKERRQ(ierr);
     300        1008 :   ierr = PetscFunctionListAdd(&SVDList,name,function);CHKERRQ(ierr);
     301        1008 :   PetscFunctionReturn(0);
     302             : }
     303             : 
     304             : /*@C
     305             :    SVDMonitorRegister - Adds SVD monitor routine.
     306             : 
     307             :    Not Collective
     308             : 
     309             :    Input Parameters:
     310             : +  name    - name of a new monitor routine
     311             : .  vtype   - a PetscViewerType for the output
     312             : .  format  - a PetscViewerFormat for the output
     313             : .  monitor - monitor routine
     314             : .  create  - creation routine, or NULL
     315             : -  destroy - destruction routine, or NULL
     316             : 
     317             :    Notes:
     318             :    SVDMonitorRegister() may be called multiple times to add several user-defined monitors.
     319             : 
     320             :    Sample usage:
     321             : .vb
     322             :    SVDMonitorRegister("my_monitor",PETSCVIEWERASCII,PETSC_VIEWER_ASCII_INFO_DETAIL,MyMonitor,NULL,NULL);
     323             : .ve
     324             : 
     325             :    Then, your monitor can be chosen with the procedural interface via
     326             : $      SVDMonitorSetFromOptions(svd,"-svd_monitor_my_monitor","my_monitor",NULL)
     327             :    or at runtime via the option
     328             : $      -svd_monitor_my_monitor
     329             : 
     330             :    Level: advanced
     331             : 
     332             : .seealso: SVDMonitorRegisterAll()
     333             : @*/
     334         864 : PetscErrorCode SVDMonitorRegister(const char name[],PetscViewerType vtype,PetscViewerFormat format,PetscErrorCode (*monitor)(SVD,PetscInt,PetscInt,PetscReal*,PetscReal*,PetscInt,PetscViewerAndFormat*),PetscErrorCode (*create)(PetscViewer,PetscViewerFormat,void*,PetscViewerAndFormat**),PetscErrorCode (*destroy)(PetscViewerAndFormat**))
     335             : {
     336         864 :   char           key[PETSC_MAX_PATH_LEN];
     337         864 :   PetscErrorCode ierr;
     338             : 
     339         864 :   PetscFunctionBegin;
     340         864 :   ierr = SVDInitializePackage();CHKERRQ(ierr);
     341         864 :   ierr = SlepcMonitorMakeKey_Internal(name,vtype,format,key);CHKERRQ(ierr);
     342         864 :   ierr = PetscFunctionListAdd(&SVDMonitorList,key,monitor);CHKERRQ(ierr);
     343         864 :   if (create)  { ierr = PetscFunctionListAdd(&SVDMonitorCreateList,key,create);CHKERRQ(ierr); }
     344         864 :   if (destroy) { ierr = PetscFunctionListAdd(&SVDMonitorDestroyList,key,destroy);CHKERRQ(ierr); }
     345         864 :   PetscFunctionReturn(0);
     346             : }
     347             : 
     348             : /*@
     349             :    SVDSetBV - Associates basis vectors objects to the singular value solver.
     350             : 
     351             :    Collective on svd
     352             : 
     353             :    Input Parameters:
     354             : +  svd - singular value solver context obtained from SVDCreate()
     355             : .  V   - the basis vectors object for right singular vectors
     356             : -  U   - the basis vectors object for left singular vectors
     357             : 
     358             :    Note:
     359             :    Use SVDGetBV() to retrieve the basis vectors contexts (for example,
     360             :    to free them at the end of the computations).
     361             : 
     362             :    Level: advanced
     363             : 
     364             : .seealso: SVDGetBV()
     365             : @*/
     366           3 : PetscErrorCode SVDSetBV(SVD svd,BV V,BV U)
     367             : {
     368           3 :   PetscErrorCode ierr;
     369             : 
     370           3 :   PetscFunctionBegin;
     371           3 :   PetscValidHeaderSpecific(svd,SVD_CLASSID,1);
     372           3 :   if (V) {
     373           3 :     PetscValidHeaderSpecific(V,BV_CLASSID,2);
     374           3 :     PetscCheckSameComm(svd,1,V,2);
     375           3 :     ierr = PetscObjectReference((PetscObject)V);CHKERRQ(ierr);
     376           3 :     ierr = BVDestroy(&svd->V);CHKERRQ(ierr);
     377           3 :     svd->V = V;
     378           3 :     ierr = PetscLogObjectParent((PetscObject)svd,(PetscObject)svd->V);CHKERRQ(ierr);
     379             :   }
     380           3 :   if (U) {
     381           3 :     PetscValidHeaderSpecific(U,BV_CLASSID,3);
     382           3 :     PetscCheckSameComm(svd,1,U,3);
     383           3 :     ierr = PetscObjectReference((PetscObject)U);CHKERRQ(ierr);
     384           3 :     ierr = BVDestroy(&svd->U);CHKERRQ(ierr);
     385           3 :     svd->U = U;
     386           3 :     ierr = PetscLogObjectParent((PetscObject)svd,(PetscObject)svd->U);CHKERRQ(ierr);
     387             :   }
     388           3 :   PetscFunctionReturn(0);
     389             : }
     390             : 
     391             : /*@
     392             :    SVDGetBV - Obtain the basis vectors objects associated to the singular
     393             :    value solver object.
     394             : 
     395             :    Not Collective
     396             : 
     397             :    Input Parameter:
     398             : .  svd - singular value solver context obtained from SVDCreate()
     399             : 
     400             :    Output Parameters:
     401             : +  V - basis vectors context for right singular vectors
     402             : -  U - basis vectors context for left singular vectors
     403             : 
     404             :    Level: advanced
     405             : 
     406             : .seealso: SVDSetBV()
     407             : @*/
     408         290 : PetscErrorCode SVDGetBV(SVD svd,BV *V,BV *U)
     409             : {
     410         290 :   PetscErrorCode ierr;
     411             : 
     412         290 :   PetscFunctionBegin;
     413         290 :   PetscValidHeaderSpecific(svd,SVD_CLASSID,1);
     414         290 :   if (V) {
     415         146 :     if (!svd->V) {
     416         146 :       ierr = BVCreate(PetscObjectComm((PetscObject)svd),&svd->V);CHKERRQ(ierr);
     417         146 :       ierr = PetscObjectIncrementTabLevel((PetscObject)svd->V,(PetscObject)svd,0);CHKERRQ(ierr);
     418         146 :       ierr = PetscLogObjectParent((PetscObject)svd,(PetscObject)svd->V);CHKERRQ(ierr);
     419         146 :       ierr = PetscObjectSetOptions((PetscObject)svd->V,((PetscObject)svd)->options);CHKERRQ(ierr);
     420             :     }
     421         146 :     *V = svd->V;
     422             :   }
     423         290 :   if (U) {
     424         144 :     if (!svd->U) {
     425         144 :       ierr = BVCreate(PetscObjectComm((PetscObject)svd),&svd->U);CHKERRQ(ierr);
     426         144 :       ierr = PetscObjectIncrementTabLevel((PetscObject)svd->U,(PetscObject)svd,0);CHKERRQ(ierr);
     427         144 :       ierr = PetscLogObjectParent((PetscObject)svd,(PetscObject)svd->U);CHKERRQ(ierr);
     428         144 :       ierr = PetscObjectSetOptions((PetscObject)svd->U,((PetscObject)svd)->options);CHKERRQ(ierr);
     429             :     }
     430         144 :     *U = svd->U;
     431             :   }
     432         290 :   PetscFunctionReturn(0);
     433             : }
     434             : 
     435             : /*@
     436             :    SVDSetDS - Associates a direct solver object to the singular value solver.
     437             : 
     438             :    Collective on svd
     439             : 
     440             :    Input Parameters:
     441             : +  svd - singular value solver context obtained from SVDCreate()
     442             : -  ds  - the direct solver object
     443             : 
     444             :    Note:
     445             :    Use SVDGetDS() to retrieve the direct solver context (for example,
     446             :    to free it at the end of the computations).
     447             : 
     448             :    Level: advanced
     449             : 
     450             : .seealso: SVDGetDS()
     451             : @*/
     452           0 : PetscErrorCode SVDSetDS(SVD svd,DS ds)
     453             : {
     454           0 :   PetscErrorCode ierr;
     455             : 
     456           0 :   PetscFunctionBegin;
     457           0 :   PetscValidHeaderSpecific(svd,SVD_CLASSID,1);
     458           0 :   PetscValidHeaderSpecific(ds,DS_CLASSID,2);
     459           0 :   PetscCheckSameComm(svd,1,ds,2);
     460           0 :   ierr = PetscObjectReference((PetscObject)ds);CHKERRQ(ierr);
     461           0 :   ierr = DSDestroy(&svd->ds);CHKERRQ(ierr);
     462           0 :   svd->ds = ds;
     463           0 :   ierr = PetscLogObjectParent((PetscObject)svd,(PetscObject)svd->ds);CHKERRQ(ierr);
     464           0 :   PetscFunctionReturn(0);
     465             : }
     466             : 
     467             : /*@
     468             :    SVDGetDS - Obtain the direct solver object associated to the singular value
     469             :    solver object.
     470             : 
     471             :    Not Collective
     472             : 
     473             :    Input Parameters:
     474             : .  svd - singular value solver context obtained from SVDCreate()
     475             : 
     476             :    Output Parameter:
     477             : .  ds - direct solver context
     478             : 
     479             :    Level: advanced
     480             : 
     481             : .seealso: SVDSetDS()
     482             : @*/
     483         143 : PetscErrorCode SVDGetDS(SVD svd,DS *ds)
     484             : {
     485         143 :   PetscErrorCode ierr;
     486             : 
     487         143 :   PetscFunctionBegin;
     488         143 :   PetscValidHeaderSpecific(svd,SVD_CLASSID,1);
     489         143 :   PetscValidPointer(ds,2);
     490         143 :   if (!svd->ds) {
     491         143 :     ierr = DSCreate(PetscObjectComm((PetscObject)svd),&svd->ds);CHKERRQ(ierr);
     492         143 :     ierr = PetscObjectIncrementTabLevel((PetscObject)svd->ds,(PetscObject)svd,0);CHKERRQ(ierr);
     493         143 :     ierr = PetscLogObjectParent((PetscObject)svd,(PetscObject)svd->ds);CHKERRQ(ierr);
     494         143 :     ierr = PetscObjectSetOptions((PetscObject)svd->ds,((PetscObject)svd)->options);CHKERRQ(ierr);
     495             :   }
     496         143 :   *ds = svd->ds;
     497         143 :   PetscFunctionReturn(0);
     498             : }
     499             : 

Generated by: LCOV version 1.14