FOSMAnalysis.cppGo to the documentation of this file.00001 /* ****************************************************************** ** 00002 ** OpenSees - Open System for Earthquake Engineering Simulation ** 00003 ** Pacific Earthquake Engineering Research Center ** 00004 ** ** 00005 ** ** 00006 ** (C) Copyright 2001, The Regents of the University of California ** 00007 ** All Rights Reserved. ** 00008 ** ** 00009 ** Commercial use of this program without express permission of the ** 00010 ** University of California, Berkeley, is strictly prohibited. See ** 00011 ** file 'COPYRIGHT' in main directory for information on usage and ** 00012 ** redistribution, and for a DISCLAIMER OF ALL WARRANTIES. ** 00013 ** ** 00014 ** Developed by: ** 00015 ** Frank McKenna (fmckenna@ce.berkeley.edu) ** 00016 ** Gregory L. Fenves (fenves@ce.berkeley.edu) ** 00017 ** Filip C. Filippou (filippou@ce.berkeley.edu) ** 00018 ** ** 00019 ** Reliability module developed by: ** 00020 ** Terje Haukaas (haukaas@ce.berkeley.edu) ** 00021 ** Armen Der Kiureghian (adk@ce.berkeley.edu) ** 00022 ** ** 00023 ** ****************************************************************** */ 00024 00025 // $Revision: 1.1 $ 00026 // $Date: 2003/10/27 23:45:41 $ 00027 // $Source: /usr/local/cvs/OpenSees/SRC/reliability/analysis/analysis/FOSMAnalysis.cpp,v $ 00028 00029 00030 // 00031 // Written by Terje Haukaas (haukaas@ce.berkeley.edu) 00032 // 00033 00034 #include <FOSMAnalysis.h> 00035 #include <ReliabilityAnalysis.h> 00036 #include <ReliabilityDomain.h> 00037 #include <RandomVariablePositioner.h> 00038 #include <GFunEvaluator.h> 00039 #include <GradGEvaluator.h> 00040 #include <Matrix.h> 00041 #include <Vector.h> 00042 #include <tcl.h> 00043 00044 #include <fstream> 00045 #include <iomanip> 00046 #include <iostream> 00047 using std::ifstream; 00048 using std::ios; 00049 using std::setw; 00050 using std::setprecision; 00051 using std::setiosflags; 00052 00053 00054 FOSMAnalysis::FOSMAnalysis(ReliabilityDomain *passedReliabilityDomain, 00055 GFunEvaluator *passedGFunEvaluator, 00056 GradGEvaluator *passedGradGEvaluator, 00057 Tcl_Interp *passedTclInterp, 00058 TCL_Char *passedFileName) 00059 :ReliabilityAnalysis() 00060 { 00061 theReliabilityDomain = passedReliabilityDomain; 00062 theGFunEvaluator = passedGFunEvaluator; 00063 theGradGEvaluator = passedGradGEvaluator; 00064 theTclInterp = passedTclInterp; 00065 fileName = new char[256]; 00066 strcpy(fileName,passedFileName); 00067 } 00068 00069 00070 FOSMAnalysis::~FOSMAnalysis() 00071 { 00072 if (fileName != 0) 00073 delete [] fileName; 00074 } 00075 00076 00077 00078 int 00079 FOSMAnalysis::analyze(void) 00080 { 00081 00082 // Alert the user that the FORM analysis has started 00083 opserr << "FOSM Analysis is running ... " << endln; 00084 00085 00086 // Initial declarations 00087 int i,j,k; 00088 00089 00090 // Open output file 00091 ofstream outputFile( fileName, ios::out ); 00092 00093 00094 // Get number of random variables and limit-state function 00095 int nrv = theReliabilityDomain->getNumberOfRandomVariables(); 00096 int numLsf = theReliabilityDomain->getNumberOfLimitStateFunctions(); 00097 00098 00099 // Get mean point 00100 RandomVariable *aRandomVariable; 00101 Vector meanVector(nrv); 00102 for (i=1; i<=nrv; i++) 00103 { 00104 aRandomVariable = theReliabilityDomain->getRandomVariablePtr(i); 00105 if (aRandomVariable == 0) { 00106 opserr << "FOSMAnalysis::analyze() -- Could not find" << endln 00107 << " random variable with tag #" << i << "." << endln; 00108 return -1; 00109 } 00110 meanVector(i-1) = aRandomVariable->getMean(); 00111 } 00112 00113 00114 // Establish vector of standard deviations 00115 Vector stdvVector(nrv); 00116 for (i=1; i<=nrv; i++) 00117 { 00118 aRandomVariable = theReliabilityDomain->getRandomVariablePtr(i); 00119 stdvVector(i-1) = aRandomVariable->getStdv(); 00120 } 00121 00122 00123 // Evaluate limit-state functions 00124 Vector meanEstimates(numLsf); 00125 int result; 00126 int lsf; 00127 result = theGFunEvaluator->runGFunAnalysis(meanVector); 00128 if (result < 0) { 00129 opserr << "FOSMAnalysis::analyze() - " << endln 00130 << " could not run analysis to evaluate limit-state function. " << endln; 00131 return -1; 00132 } 00133 for (lsf=1; lsf<=numLsf; lsf++ ) { 00134 theReliabilityDomain->setTagOfActiveLimitStateFunction(lsf); 00135 result = theGFunEvaluator->evaluateG(meanVector); 00136 if (result < 0) { 00137 opserr << "FOSMAnalysis::analyze() - " << endln 00138 << " could not tokenize limit-state function. " << endln; 00139 return -1; 00140 } 00141 meanEstimates(lsf-1) = theGFunEvaluator->getG(); 00142 } 00143 00144 00145 // Evaluate the gradients of limit-state functions 00146 Matrix matrixOfGradientVectors(nrv,numLsf); 00147 result = theGradGEvaluator->computeAllGradG(meanEstimates,meanVector); 00148 if (result < 0) { 00149 opserr << "FOSMAnalysis::analyze() -- could not" << endln 00150 << " compute gradients of the limit-state function. " << endln; 00151 return -1; 00152 } 00153 matrixOfGradientVectors = theGradGEvaluator->getAllGradG(); 00154 00155 00156 // Establish covariance matrix 00157 Matrix covMatrix(nrv,nrv); 00158 for (i=1; i<=nrv; i++) { 00159 covMatrix(i-1,i-1) = stdvVector(i-1)*stdvVector(i-1); 00160 } 00161 int ncorr = theReliabilityDomain->getNumberOfCorrelationCoefficients(); 00162 CorrelationCoefficient *theCorrelationCoefficient; 00163 double covariance, correlation; 00164 int rv1, rv2; 00165 for (i=1 ; i<=ncorr ; i++) { 00166 theCorrelationCoefficient = theReliabilityDomain->getCorrelationCoefficientPtr(i); 00167 correlation = theCorrelationCoefficient->getCorrelation(); 00168 rv1 = theCorrelationCoefficient->getRv1(); 00169 rv2 = theCorrelationCoefficient->getRv2(); 00170 covariance = correlation*stdvVector(rv1-1)*stdvVector(rv2-1); 00171 covMatrix(rv1-1,rv2-1) = covariance; 00172 covMatrix(rv2-1,rv1-1) = covariance; 00173 } 00174 00175 00176 // Post-processing loop over limit-state functions 00177 Vector responseStdv(numLsf); 00178 Vector gradient(nrv); 00179 double responseVariance; 00180 for (lsf=1; lsf<=numLsf; lsf++ ) { 00181 00182 00183 // Set tag of active limit-state function 00184 theReliabilityDomain->setTagOfActiveLimitStateFunction(lsf); 00185 00186 00187 // Extract relevant gradient 00188 for (i=0; i<nrv; i++) { 00189 gradient(i) = matrixOfGradientVectors(i,lsf-1); 00190 } 00191 00192 00193 // Estimate of standard deviation of response 00194 responseVariance = (covMatrix^gradient)^gradient; 00195 if (responseVariance <= 0.0) { 00196 opserr << "ERROR: Response variance of limit-state function number "<< lsf << endln 00197 << " is zero! " << endln; 00198 } 00199 else { 00200 responseStdv(lsf-1) = sqrt(responseVariance); 00201 } 00202 00203 // Compute importance measure (dgdx*stdv) 00204 Vector importance(nrv); 00205 for (i=1 ; i<=nrv ; i++) { 00206 aRandomVariable = theReliabilityDomain->getRandomVariablePtr(i); 00207 importance(i-1) = gradient(i-1) * aRandomVariable->getStdv(); 00208 } 00209 double imptNorm = importance.Norm(); 00210 for (i=1 ; i<=nrv ; i++) { 00211 importance(i-1) = importance(i-1)/imptNorm; 00212 } 00213 00214 00215 // Print FOSM results to the output file 00216 outputFile << "#######################################################################" << endln; 00217 outputFile << "# FOSM ANALYSIS RESULTS, LIMIT-STATE FUNCTION NUMBER " 00218 <<setiosflags(ios::left)<<setprecision(1)<<setw(4)<<lsf <<" #" << endln; 00219 outputFile << "# #" << endln; 00220 00221 outputFile << "# Estimated mean: .................................... " 00222 <<setiosflags(ios::left)<<setprecision(5)<<setw(12)<<meanEstimates(lsf-1) 00223 << " #" << endln; 00224 outputFile << "# Estimated standard deviation: ...................... " 00225 <<setiosflags(ios::left)<<setprecision(5)<<setw(12)<<responseStdv(lsf-1) 00226 << " #" << endln; 00227 outputFile << "# #" << endln; 00228 outputFile << "# R.v. Importance measure (dgdx*stdv) #" << endln; 00229 outputFile.setf( ios::scientific, ios::floatfield ); 00230 for (int i=0; i<nrv; i++) { 00231 outputFile << "# " <<setw(3)<<(i+1)<<" "; 00232 outputFile << ""; 00233 if (importance(i)<0.0) { 00234 outputFile << "-"; 00235 } 00236 else { 00237 outputFile << " "; 00238 } 00239 outputFile <<setprecision(3)<<setw(11)<<fabs(importance(i)); 00240 outputFile << " #" << endln; 00241 } 00242 outputFile << "# #" << endln; 00243 outputFile << "#######################################################################" << endln << endln << endln; 00244 outputFile.flush(); 00245 00246 } 00247 00248 00249 // Estimation of response covariance matrix 00250 Matrix responseCovMatrix(numLsf,numLsf); 00251 double responseCovariance; 00252 Vector gradientVector1(nrv), gradientVector2(nrv); 00253 for (i=1; i<=numLsf; i++) { 00254 for (k=0; k<nrv; k++) { 00255 gradientVector1(k) = matrixOfGradientVectors(k,i-1); 00256 } 00257 for (j=i+1; j<=numLsf; j++) { 00258 for (k=0; k<nrv; k++) { 00259 gradientVector2(k) = matrixOfGradientVectors(k,j-1); 00260 } 00261 responseCovariance = (covMatrix^gradientVector1)^gradientVector2; 00262 responseCovMatrix(i-1,j-1) = responseCovariance; 00263 } 00264 } 00265 for (i=1; i<=numLsf; i++) { 00266 for (j=1; j<i; j++) { 00267 responseCovMatrix(i-1,j-1) = responseCovMatrix(j-1,i-1); 00268 } 00269 } 00270 00271 00272 // Corresponding correlation matrix 00273 Matrix correlationMatrix(numLsf,numLsf); 00274 for (i=1; i<=numLsf; i++) { 00275 for (j=i+1; j<=numLsf; j++) { 00276 correlationMatrix(i-1,j-1) = responseCovMatrix(i-1,j-1)/(responseStdv(i-1)*responseStdv(j-1)); 00277 } 00278 } 00279 for (i=1; i<=numLsf; i++) { 00280 for (j=1; j<i; j++) { 00281 correlationMatrix(i-1,j-1) = correlationMatrix(j-1,i-1); 00282 } 00283 } 00284 00285 00286 // Print correlation results 00287 outputFile << "#######################################################################" << endln; 00288 outputFile << "# RESPONSE CORRELATION COEFFICIENTS #" << endln; 00289 outputFile << "# #" << endln; 00290 if (numLsf <=1) { 00291 outputFile << "# Only one limit-state function! #" << endln; 00292 } 00293 else { 00294 outputFile << "# gFun gFun Correlation #" << endln; 00295 outputFile.setf(ios::fixed, ios::floatfield); 00296 for (i=0; i<numLsf; i++) { 00297 for (j=i+1; j<numLsf; j++) { 00298 // outputFile.setf(ios::fixed, ios::floatfield); 00299 outputFile << "# " <<setw(3)<<(i+1)<<" "<<setw(3)<<(j+1)<<" "; 00300 if (correlationMatrix(i,j)<0.0) { outputFile << "-"; } 00301 else { outputFile << " "; } 00302 // outputFile.setf(ios::scientific, ios::floatfield); 00303 outputFile <<setprecision(7)<<setw(11)<<fabs(correlationMatrix(i,j)); 00304 outputFile << " #" << endln; 00305 } 00306 } 00307 } 00308 outputFile << "# #" << endln; 00309 outputFile << "#######################################################################" << endln << endln << endln; 00310 00311 00312 00313 // Print summary of results to screen (more here!!!) 00314 opserr << "FOSMAnalysis completed." << endln; 00315 00316 00317 return 0; 00318 } 00319
Generated on Mon Oct 23 15:05:25 2006 for OpenSees by |
opensees-support @ berkeley.edu
Supported by the National Science Foundation
©2006, UC Regents