CurvaturesBySearchAlgorithm.cpp

Go 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.4 $
00026 // $Date: 2003/03/04 00:38:55 $
00027 // $Source: /usr/local/cvs/OpenSees/SRC/reliability/analysis/curvature/CurvaturesBySearchAlgorithm.cpp,v $
00028 
00029 
00030 //
00031 // Written by Terje Haukaas (haukaas@ce.berkeley.edu) 
00032 //
00033 
00034 #include <CurvaturesBySearchAlgorithm.h>
00035 #include <FindCurvatures.h>
00036 #include <LimitStateFunction.h>
00037 #include <FindDesignPointAlgorithm.h>
00038 #include <RandomVariable.h>
00039 #include <Vector.h>
00040 
00041 
00042 CurvaturesBySearchAlgorithm::CurvaturesBySearchAlgorithm(int passedNumberOfCurvatures,
00043                                                                         FindDesignPointAlgorithm *passedFindDesignPointAlgorithm)
00044 :FindCurvatures(), curvatures(passedNumberOfCurvatures)
00045 {
00046         numberOfCurvatures = passedNumberOfCurvatures;
00047         theFindDesignPointAlgorithm = passedFindDesignPointAlgorithm;
00048 }
00049 
00050 CurvaturesBySearchAlgorithm::~CurvaturesBySearchAlgorithm()
00051 {
00052 }
00053 
00054 
00055 int
00056 CurvaturesBySearchAlgorithm::computeCurvatures(ReliabilityDomain *theReliabilityDomain)
00057 {
00058 
00059         // Initial declarations
00060         Vector last_u;
00061         Vector secondLast_u;
00062         Vector lastAlpha;
00063         Vector secondLastAlpha;
00064         Vector lastDirection;
00065         Vector uLastMinus_u = last_u - secondLast_u;
00066         double signumProduct;
00067         double alphaProduct;
00068         double sumSquared;
00069         double norm_uLastMinus_u;
00070 
00071         // "Download" limit-state function from reliability domain
00072         int lsf = theReliabilityDomain->getTagOfActiveLimitStateFunction();
00073         LimitStateFunction *theLimitStateFunction = 
00074                 theReliabilityDomain->getLimitStateFunctionPtr(lsf);
00075 
00076         // The design point in the original space
00077         Vector x_star = theLimitStateFunction->designPoint_x_inOriginalSpace;
00078 
00079         // Number of random variables
00080         int nrv = x_star.Size();
00081 
00082         // Declare vector to store all principal axes
00083         Vector principalAxes( nrv * numberOfCurvatures );
00084 
00085         // Start point of new searches: Perturb x^star by 10% of each standard deviation
00086         Vector startPoint(nrv);
00087         RandomVariable *aRandomVariable;
00088         int numberOfRandomVariables = theReliabilityDomain->getNumberOfRandomVariables();
00089         int i;
00090         for (i=0; i<nrv; i++) {
00091                 aRandomVariable = theReliabilityDomain->getRandomVariablePtr(i+1);
00092                 startPoint(i) = aRandomVariable->getStartValue();
00093         }
00094 
00095         // Compute curvatures by repeated searches
00096         for (i=0; i<numberOfCurvatures; i++) {
00097 
00098                 // Get hold of 'u', 'alpha' and search direction at the two last steps
00099                 last_u = theLimitStateFunction->designPoint_u_inStdNormalSpace;
00100                 secondLast_u = theLimitStateFunction->secondLast_u;
00101                 lastAlpha = theLimitStateFunction->normalizedNegativeGradientVectorAlpha;
00102                 secondLastAlpha = theLimitStateFunction->secondLastAlpha;
00103                 lastDirection = theLimitStateFunction->lastSearchDirection;
00104 
00105                 // Compute curvature according to Der Kiureghian & De Stefano (1992), Eq.26:
00106 
00107                 // Initial computations
00108                 uLastMinus_u = last_u - secondLast_u;
00109                 signumProduct = secondLastAlpha ^ uLastMinus_u;
00110                 alphaProduct = secondLastAlpha ^ lastAlpha;
00111                 sumSquared = 0.0;
00112 
00113                 // Compute norm of the difference vector
00114                 for ( int k=0; k<last_u.Size(); k++ ) {
00115                         sumSquared += uLastMinus_u(k)*uLastMinus_u(k);
00116                 }
00117                 norm_uLastMinus_u = sqrt(sumSquared);
00118 
00119                 // Check sign and compute curvature
00120                 if (fabs(signumProduct)==(signumProduct)) {
00121                         curvatures(i) = acos(alphaProduct) / norm_uLastMinus_u;
00122                 }
00123                 else {
00124                         curvatures(i) = -acos(alphaProduct) / norm_uLastMinus_u;
00125                 }
00126 
00127                 // Run a new search in the subspace orthogonal to previous principal directions
00128                 if ( i!=(numberOfCurvatures-1) ) {
00129 
00130                         // Store all previous principal axes in a vector
00131                         for (int j=0; j<nrv; j++ ) {
00132                                 principalAxes((i*nrv)+j) = lastDirection(j);
00133                         }
00134 
00135                         // To be completed...
00136                 }
00137         }
00138 
00139         return 0;
00140 }
00141 
00142 
00143 
00144 Vector
00145 CurvaturesBySearchAlgorithm::getCurvatures()
00146 {
00147         return curvatures;
00148 }
00149 
00150 
00151 

Generated on Mon Oct 23 15:05:25 2006 for OpenSees by doxygen 1.5.0