F08 ID: 844150
Download Pdf The PPT/PDF document "NAGLibraryRoutineDocumentF08KPFZGESVDbef..." is the property of its rightful owner. Permission is granted to download and print the materials on this web site for personal, non-commercial use only, and to display it on your personal computer provided you do not modify the materials and that you retain all copyright notices contained in the materials. By downloading content from our website, you accept the terms of this agreement.
1 NAGLibraryRoutineDocumentF08KPF(ZGESVD)b
NAGLibraryRoutineDocumentF08KPF(ZGESVD)beforeusingthisroutine,pleasereadtheUsersNoteforyourimplementationtochecktheinterpretationofbolditalicisedtermsandotherimplementation-dependentdetails.1PurposeF08KPF(ZGESVD)computesthesingularvaluedecomposition(SVD)ofacomplex F08 Least-squaresandEigenvalueProblems(LAPACK)Mark24F08KPF.1 Nocolumnsof(noleftsingularvectors)arecomputed.:JOBUA,S,OorN.JOBVT CHARACTER(1)Onentry:speciesoptionsforcomputingallorpartofthematrixrowsofarereturnedinthearrayVT.Therstminm;nrowsof(therightsingularvectors)arereturnedinthearrayVT.Therstminm;nrowsof(therightsingularvectors)areoverwrittenonthearrayA.Norowsof(norightsingularvectors)arecomputed.A,S,OorN;JOBVTandJOBUcannotbothbeO.M INTEGEROnentry,thenumberofrowsofthematrixN INTEGEROnentry,thenumberofcolumnsofthematrix COMPLEX(KIND=nag_wp)array:theseconddimensionofthearrayAmustbeatleastmax1Onentry:theOnexit:ifJOBUO,Aisoverwrittenwiththerstminm;ncolumnsof(theleftsingularvectors,storedcolumn-wise).IfJOBVTO,Aisoverwrittenwiththerstminm;nrowsof(therightsingularvectors,storedrow-wise).IfJOBUOandJOBVTO,thecontentsofAaredestroyed.LDA INTEGEROnentry:therstdimensionofthearrayAasdeclaredinthe(sub)programfromwhichF08KPF(ZGESVD)iscalled.:LDAmax1
2 REAL(KIND=nag_wp)array:thedimensionofth
REAL(KIND=nag_wp)array:thedimensionofthearraySmustbeatleastmax1minMOnexit:thesingularvaluesof,sortedsothatS NAGLibraryManualF08KPF.2Mark24 COMPLEX(KIND=nag_wp)array:theseconddimensionofthearrayUmustbeatleastmax1ifJOBUmax1minMifJOBUS,andatleast1otherwise.Onexit:ifJOBUA,UcontainstheunitarymatrixIfJOBUS,Ucontainstherstminm;ncolumnsof(theleftsingularvectors,storedIfJOBUNorO,Uisnotreferenced.LDU INTEGEROnentry:therstdimensionofthearrayUasdeclaredinthe(sub)programfromwhichF08KPF(ZGESVD)iscalled.ifJOBUAorS,LDUmax1otherwiseLDULDVT, COMPLEX(KIND=nag_wp)array:theseconddimensionofthearrayVTmustbeatleastmax1ifJOBVTAorS,andatleast1otherwise.Onexit:ifJOBVTA,VTcontainstheunitarymatrixIfJOBVTS,VTcontainstherstminm;nrowsof(therightsingularvectors,storedrow-IfJOBVTNorO,VTisnotreferenced.11:LDVT INTEGEROnentry:therstdimensionofthearrayVTasdeclaredinthe(sub)programfromwhichF08KPF(ZGESVD)iscalled.ifJOBVTA,LDVTmax1ifJOBVTS,LDVTmax1minMotherwiseLDVTmax1ðÞÞ COMPLEX(KIND=nag_wp)arrayWorkspaceOnexit:ifINFO0,therealpartofWORKcontainstheminimumvalueofLWORKrequiredforoptimalperformance.LWORK INTEGEROnentry:thedimensionofthearrayWORKasdeclaredinthe(sub)programfromwhichF08KPF(ZGESVD)iscalled.IfLWORK1,aworkspacequeryisassumed;th
3 eroutineonlycalculatestheoptimalsizeofth
eroutineonlycalculatestheoptimalsizeoftheWORKarray,returnsthisvalueastherstentryoftheWORKarray,andnoerrormessagerelatedtoLWORKisissued.Suggestedvalue:foroptimalperformance,LWORKshouldgenerallybelarger.ConsiderincreasingLWORKbyatleastminM,whereistheoptimalblocksize:LWORKmax1minMðÞþmaxM F08 Least-squaresandEigenvalueProblems(LAPACK)Mark24F08KPF.3 RWORK REAL(KIND=nag_wp)arrayWorkspace:thedimensionofthearrayRWORKmustbeatleastmax1minMOnexit:ifINFO0,RWORKminMðÞcontainstheunconvergedsuperdiagonalelementsofanupperbidiagonalmatrixwhosediagonalisin(notnecessarilysorted).UBV,soithasthesamesingularvaluesas,andsingularvectorsrelatedbyINFO INTEGEROnexit:INFO0unlesstheroutinedetectsanerror(seeSection6).6ErrorIndicatorsandWarningsErrorsorwarningsdetectedbytheroutine:IfINFO,argumenthadanillegalvalue.Anexplanatorymessageisoutput,andexecutionoftheprogramisterminated.IfF08KPF(ZGESVD)didnotconverge,INFOspecieshowmanysuperdiagonalsofanintermediatebidiagonalformdidnotconvergetozero.SeethedescriptionofRWORKabovefor7AccuracyThecomputedsingularvaluedecompositionisnearlytheexactsingularvaluedecompositionforanearby,whereisthemachineprecision.Inaddition,thecomputedsingularvectorsarenearlyorthogonaltoworkingprecision.SeeSection4.
4 9ofAndersonetal.(1999)forfurtherdetails.
9ofAndersonetal.(1999)forfurtherdetails.8FurtherCommentsThetotalnumberofoatingpointoperationsisapproximatelyproportionaltoThesingularvaluesarereturnedindescendingorder.TherealanalogueofthisroutineisF08KBF(DGESVD).9ExampleThisexamplendsthesingularvaluesandleftandrightsingularvectorsofthe6by4matrixtogetherwithapproximateerrorboundsforthecomputedsingularvaluesandvectors.TheexampleprogramforF08KRF(ZGESDD)illustratesndingasingularvaluedecompositionfortheNAGLibraryManualF08KPF.4Mark24 9.1ProgramTextProgramf08kpfe!F08KPFExampleProgramText!Mark24Release.NAGCopyright2012.!..UseStatements..Usenag_library,Only:nag_wp,zgesvd!..ImplicitNoneStatement..ImplicitNone!..Parameters..Integer,Parameter::nb=64,nin=5,nout=6,&prerr=0!..LocalScalars..Integer::i,info,lda,ldu,ldvt,&lwork,m,n!..LocalArrays..Complex(Kind=nag_wp),Allocatable::a(:,:),a_copy(:,:),b(:),&u(:,:),vt(:,:),work(:)Complex(Kind=nag_wp)::dummy(1,1)Real(Kind=nag_wp),Allocatable::rwork(:),s(:)!..IntrinsicProcedures..Intrinsic::max,min,nint,real!..ExecutableStatements..ContinueWrite(nout,)F08KPFExampleProgramResultsWrite(nout,!SkipheadingindatafileRead(nin,Read(nin,)m,nlda=mldu=mldvt=nAllocate(a(lda,n),a_copy(m,n),s(n),u(ldu,m),vt(ldvt,n),b(m),rwork(5!Readth
5 embynmatrixAfromdatafileRead(nin,)(a(i,1
embynmatrixAfromdatafileRead(nin,)(a(i,1:n),i=1,m)!ReadtherighthandsideofthelinearsystemRead(nin,)b(1:m)a_copy(1:m,1:n)=a(1:m,1:n)!Useroutineworkspacequerytogetoptimalworkspace.lwork=-1!TheNAGnameequivalentofdgesvdisf08kpfCallzgesvd(A,S,m,n,a,lda,s,u,ldu,vt,ldvt,dummy,lwork,rwork,info)!Makesurethatthereisenoughworkspaceforblocksizenb.lwork=max(m+3(m+n),nint(real(dummy(1,1))))Allocate(work(lwork))!Computethesingularvaluesandleftandrightsingularvectors!ofA.!TheNAGnameequivalentofdgesvdisf08kpfCallzgesvd(A,S,m,n,a,lda,s,u,ldu,vt,ldvt,work,lwork,rwork,info)If(info/=0)ThenWrite(nout,99999)FailureinF08KPF/ZGESVD.INFO=,info99999Format(1X,A,I4)GoTo100EndIf!PrintthesignificantsingularvaluesofAWrite(nout,)SingularvaluesofA:Write(nout,99998)s(1:min(m,n))99998Format(1X,4(3X,F11.4))F08 Least-squaresandEigenvalueProblems(LAPACK)Mark24F08KPF.5 If�(prerr0)ThenCallcompute_error_bounds(m,n,s)EndIfIf�(mn)Then!ComputeVbtogetleast-squaressolution.Callcompute_least_squares(m,n,a_copy,m,u,ldu,vt,ldvt,s,b)EndIf100ContinueContainsSubroutinecompute_least_squares(m,n,a,lda,u,ldu,vt,ldvt,s,b)!..UseStatements..Usenag_library,Only:dznrm2,zgemv!..ImplicitNoneStatement..ImplicitNone!..ScalarArguments..Integer
6 ,Intent(In)::lda,ldu,ldvt,m,n!..ArrayArg
,Intent(In)::lda,ldu,ldvt,m,n!..ArrayArguments..Complex(Kind=nag_wp),Intent(In)::a(lda,n),u(ldu,m),vt(ldvt,n)Complex(Kind=nag_wp),Intent(Inout)::b(m)Real(Kind=nag_wp),Intent(In)::s(n)!..LocalScalars..Complex(Kind=nag_wp)::alpha,betaReal(Kind=nag_wp)::norm!..LocalArrays..Complex(Kind=nag_wp),Allocatable::x(:),y(:)!..IntrinsicProcedures..Intrinsic::allocated,cmplx!..ExecutableStatements..Allocate(x(n),y(n))!ComputeVbtogetleast-squaressolution.!y=U^Tb!TheNAGnameequivelentofzgemvisf06safalpha=cmplx(1.0_nag_wp,0.0_nag_wp,kind=nag_wp)beta=cmplx(0.0_nag_wp,0.0_nag_wp,kind=nag_wp)Callzgemv(C,m,n,alpha,u,ldu,b,1,beta,y,1)y(1:n)=y(1:n)/s(1:n)!x=VyCallzgemv(C,n,n,alpha,vt,ldvt,y,1,beta,x,1)Write(nout,Write(nout,)Leastsquaressolution:Write(nout,99999)x(1:n)!Findnormofresidual||b-Ax||.alpha=cmplx(-1.0_nag_wp,0.0_nag_wp,kind=nag_wp)beta=cmplx(1._nag_wp,0.0_nag_wp,kind=nag_wp)Callzgemv(N,m,n,alpha,a,lda,x,1,beta,b,1)norm=dznrm2(m,b,1)Write(nout,Write(nout,)NormofResidual:Write(nout,99998)normIf(allocated(x))Deallocate(x)If(allocated(y))Deallocate(y)99999Format(4X,(,F8.4,,,F8.4,))99998Format(4X,F11.4)EndSubroutinecompute_least_squaresNAGLibraryManualF08KPF.6Mark24 Subroutinecompute_error_bounds(m,n,s)!E
7 rrorestimatesforsingularvaluesandvectors
rrorestimatesforsingularvaluesandvectorsiscomputed!andprintedhere.!..UseStatements..Usenag_library,Only:ddisna,nag_wp,x02ajf!..ImplicitNoneStatement..ImplicitNone!..ScalarArguments..Integer,Intent(In)::m,n!..ArrayArguments..Real(Kind=nag_wp),Intent(In)::s(n)!..LocalScalars..Real(Kind=nag_wp)::eps,serrbdInteger::i,info!..LocalArrays..Real(Kind=nag_wp),Allocatable::rcondu(:),rcondv(:),&uerrbd(:),verrbd(:)!..ExecutableStatements..Allocate(rcondu(n),rcondv(n),uerrbd(n),verrbd(n))!Getthemachineprecision,EPSandcomputetheapproximate!errorboundforthecomputedsingularvalues.Notethatfor!the2-norm,S(1)=norm(A)eps=x02ajf()serrbd=eps!CallDDISNA(F08FLF)toestimatereciprocalcondition!numbersforthesingularvectorsCallddisna(Left,m,n,s,rcondu,info)Callddisna(Right,m,n,s,rcondv,info)!ComputetheerrorestimatesforthesingularvectorsDoi=1,nuerrbd(i)=serrbd/rcondu(i)verrbd(i)=serrbd/rcondv(i)EndDo!Printtheapproximateerrorboundsforthesingularvalues!andvectorsWrite(nout,Write(nout,)ErrorestimateforthesingularvaluesWrite(nout,99999)serrbdWrite(nout,Write(nout,)ErrorestimatesfortheleftsingularvectorsWrite(nout,99999)uerrbd(1:n)Write(nout,Write(nout,)ErrorestimatesfortherightsingularvectorsWrite(nout,99999)verrbd(1:n)99999Fo
8 rmat(4X,1P,6E11.1)EndSubroutinecompute_e
rmat(4X,1P,6E11.1)EndSubroutinecompute_error_boundsEndProgramf08kpfe9.2ProgramDataF08KPFExampleProgramData64:mandn(0.96,-0.81)(-0.03,0.96)(-0.91,2.06)(-0.05,0.41)F08 Least-squaresandEigenvalueProblems(LAPACK)Mark24F08KPF.7 (-0.98,1.98)(-1.20,0.19)(-0.66,0.42)(-0.81,0.56)(0.62,-0.46)(1.01,0.02)(0.63,-0.17)(-1.11,0.60)(-0.37,0.38)(0.19,-0.54)(-0.98,-0.36)(0.22,-0.20)(0.83,0.51)(0.20,0.01)(-0.17,-0.46)(1.47,1.59)(1.08,-0.28)(0.20,-0.12)(-0.07,1.23)(0.26,0.26):MatrixA(1:m,1:n)(1.00,0.00)(1.00,0.00)(1.00,0.00)(1.00,0.00)(1.00,0.00)(1.00,0.00):RHSb(1:n)9.3ProgramResultsF08KPFExampleProgramResultsSingularvaluesofA:3.99943.00031.99440.9995Leastsquaressolution:(-0.0719,-0.2761)(0.6796,0.4326)(-0.3832,-0.5447)(0.0096,-0.3489)NormofResidual:___________________________________________________________________________________________________________________________________________________________________________________________________________________________________________________________________________________________________________________________________________________________________________________________________________________________________________________________________NAGLibraryManualF08KPF.8(last)M