Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

InverseFORM: Fix code with 1.23 #42

Merged
merged 1 commit into from
May 2, 2024
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
63 changes: 32 additions & 31 deletions lib/src/InverseFORM.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -228,43 +228,44 @@ Function InverseFORM::getG(const Scalar p)
newFunction.setParameter(params);
RandomVector antecedent(event_.getImplementation()->getAntecedent().getImplementation()->clone());
const Distribution distribution(antecedent.getDistribution());
if(distribution.getImplementation()->getClassName() == "ComposedDistribution")
#if OPENTURNS_VERSION >= 102300
const JointDistribution * p_joint = dynamic_cast<JointDistribution *>(distribution.getImplementation().get());
#else
const ComposedDistribution * p_joint = dynamic_cast<ComposedDistribution *>(distribution.getImplementation().get());
#endif
if (p_joint)
{
const ComposedDistribution * p_composed = dynamic_cast<ComposedDistribution *>(distribution.getImplementation().get());
if (p_composed)
ComposedDistribution::DistributionCollection distributionCollection(p_joint->getDistributionCollection());
for (UnsignedInteger i = 0; i < distributionCollection.getSize(); ++ i)
{
ComposedDistribution::DistributionCollection distributionCollection(p_composed->getDistributionCollection());
for (UnsignedInteger i = 0; i < distributionCollection.getSize(); ++ i)
if (distributionCollection[i].getImplementation()->getClassName() == "ConditionalDistribution")
{
if (distributionCollection[i].getImplementation()->getClassName() == "ConditionalDistribution")
DistributionImplementation::PointWithDescriptionCollection parametersCollection(distributionCollection[i].getParametersCollection());
for (UnsignedInteger j = 0; j < parametersCollection.getSize(); ++ j)
{
DistributionImplementation::PointWithDescriptionCollection parametersCollection(distributionCollection[i].getParametersCollection());
for (UnsignedInteger j = 0; j < parametersCollection.getSize(); ++ j)
{
const String marginalName(parametersCollection[j].getName());
const size_t pos = marginalName.find(parameterName_);
if (pos != String::npos) {
const String endName = marginalName.substr(pos, marginalName.size());
if (endName == parameterName_)
{
parametersCollection[j] = params;
}
const String marginalName(parametersCollection[j].getName());
const size_t pos = marginalName.find(parameterName_);
if (pos != String::npos) {
const String endName = marginalName.substr(pos, marginalName.size());
if (endName == parameterName_)
{
parametersCollection[j] = params;
}
}
const ConditionalDistribution * p_conditional = dynamic_cast<ConditionalDistribution *>(distributionCollection[i].getImplementation().get());
if (p_conditional)
{
Distribution conditioning(p_conditional->getConditioningDistribution());
conditioning.setParametersCollection(parametersCollection);
ConditionalDistribution newConditional(p_conditional->getConditionedDistribution(), conditioning);
distributionCollection[i] = newConditional;
ComposedDistribution newDistribution(distributionCollection);
antecedent = RandomVector(newDistribution);
} // if p_conditional
} // if conditional
} // i
} // if p_composed
} // if composed
}
const ConditionalDistribution * p_conditional = dynamic_cast<ConditionalDistribution *>(distributionCollection[i].getImplementation().get());
if (p_conditional)
{
Distribution conditioning(p_conditional->getConditioningDistribution());
conditioning.setParametersCollection(parametersCollection);
ConditionalDistribution newConditional(p_conditional->getConditionedDistribution(), conditioning);
distributionCollection[i] = newConditional;
ComposedDistribution newDistribution(distributionCollection);
antecedent = RandomVector(newDistribution);
} // if p_conditional
} // if conditional
} // i
} // if p_joint

const CompositeRandomVector composite(newFunction, antecedent);
const ThresholdEvent event(composite, op, threshold);
Expand Down
1 change: 1 addition & 0 deletions lib/test/CMakeLists.txt
Original file line number Diff line number Diff line change
Expand Up @@ -48,6 +48,7 @@ ot_check_test (SequentialMonteCarloRobustAlgorithm_std)
endif()
ot_check_test (SubsetInverseSampling_std)
ot_check_test (InverseFORM_std IGNOREOUT)
ot_check_test (InverseFORM_sphere IGNOREOUT)

add_custom_target ( cppcheck COMMAND ${CMAKE_CTEST_COMMAND} -R "^cppcheck_"
DEPENDS ${CHECK_TO_BE_RUN}
Expand Down
93 changes: 93 additions & 0 deletions lib/test/t_InverseFORM_sphere.cxx
Original file line number Diff line number Diff line change
@@ -0,0 +1,93 @@
#include "openturns/OT.hxx"
#include "openturns/OTtestcode.hxx"
#include "otrobopt/OTRobOpt.hxx"

using namespace OT;
using namespace OT::Test;
using namespace OTROBOPT;


int main()
{
TESTPREAMBLE;
Log::Show(Log::INFO);
OStream fullprint(std::cout);

try
{
// sphere under pressure
const UnsignedInteger dim = 3;
const SymbolicFunction function(Description({"R", "e", "mulog_e", "p"}), Description({"700.0-p*R/(2.*e)"}));

const Scalar L0 = -4.715;
ParametricFunction parametric(function, Indices({2}), Point({L0}));

Dirac mulog_eDist(L0);
mulog_eDist.setDescription(Description(1, "mulog_e"));
ComposedDistribution::DistributionCollection eColl; eColl.add(mulog_eDist); eColl.add(Dirac(0.1)); eColl.add(Dirac(0.));
ComposedDistribution eParams(eColl);

ComposedDistribution::DistributionCollection coll;
coll.add(Beta(0.117284, 0.117284, 2.9, 3.1));//R
ConditionalDistribution eDist(LogNormal(L0, 0.1, 0.), eParams);

coll.add(eDist);//e
coll.add(WeibullMin(3.16471, 9.21097, 0.0));//p
ComposedDistribution myDistribution(coll);

Point median(dim);
for(UnsignedInteger i = 0; i < dim; ++ i)
median[i] = myDistribution.getMarginal(i).computeQuantile(0.5)[0];

/* We create a 'usual' RandomVector from the Distribution */
RandomVector vect(myDistribution);

/* We create a composite random vector */
CompositeRandomVector output(parametric, vect);

/* We create an Event from this RandomVector */
ThresholdEvent event(output, Less(), 0.);

fullprint << "event=" << event << std::endl;

// create the algorithm
InverseFORM algo(event, "mulog_e", median);

// print result
algo.run();
const InverseFORMResult result(algo.getResult());
fullprint << "Inverse FORM result=" << result << std::endl;
assert_almost_equal(result.getParameter(), Point({-4.69056}), 0.0, 0.03);

// FORM must yield the same probability on the limit state with parameter set to the optimum
eColl[0] = Dirac(result.getParameter()[0]);
eParams = ComposedDistribution(eColl);
eDist = ConditionalDistribution(LogNormal(result.getParameter()[0], 0.1, 0.0), eParams);
coll[1] = eDist;
myDistribution = ComposedDistribution(coll);
vect = RandomVector(myDistribution);
parametric.setParameter(result.getParameter());
output = CompositeRandomVector(parametric, vect);
event = ThresholdEvent(output, Less(), 0.);

Cobyla cobyla;
cobyla.setMaximumIterationNumber(100);
cobyla.setMaximumAbsoluteError(1.0e-10);
cobyla.setMaximumRelativeError(1.0e-10);
cobyla.setMaximumResidualError(1.0e-10);
cobyla.setMaximumConstraintError(1.0e-10);

FORM form(cobyla, event, median);
form.run();
const FORMResult resultFORM(form.getResult());
fullprint << "result FORM="<<resultFORM << std::endl;
assert_almost_equal(resultFORM.getHasoferReliabilityIndex(), algo.getTargetBeta(), 0.0, 0.1);
}
catch (const TestFailed & ex)
{
std::cerr << ex << std::endl;
return ExitCode::Error;
}
return ExitCode::Success;
}

81 changes: 4 additions & 77 deletions lib/test/t_InverseFORM_std.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -10,16 +10,15 @@ using namespace OTROBOPT;
int main()
{
TESTPREAMBLE;

try {
Log::Show(Log::INFO);

OStream fullprint(std::cout);

try
{
// beam model
UnsignedInteger dim = 4;
const UnsignedInteger dim = 4;

SymbolicFunction function(Description({"E", "F", "L", "b", "h"}), Description({"F*L^3/(48.*E*b*h^3/12.)"}));
const SymbolicFunction function(Description({"E", "F", "L", "b", "h"}), Description({"F*L^3/(48.*E*b*h^3/12.)"}));
ParametricFunction parametric(function, Indices({2}), Point({5.0}));

ComposedDistribution::DistributionCollection coll;
Expand All @@ -30,7 +29,6 @@ int main()

const ComposedDistribution myDistribution(coll);


Point median(dim);
for(UnsignedInteger i = 0; i < dim; ++ i)
median[i] = myDistribution.getMarginal(i).computeQuantile(0.5)[0];
Expand Down Expand Up @@ -73,77 +71,6 @@ int main()
fullprint << "FORM result="<< resultFORM << std::endl;
assert_almost_equal(resultFORM.getHasoferReliabilityIndex(), algo.getTargetBeta());
}

{
// sphere under pressure
UnsignedInteger dim = 3;
SymbolicFunction function(Description({"R", "e", "mulog_e", "p"}), Description({"700.0-p*R/(2.*e)"}));

const Scalar L0 = -4.715;
ParametricFunction parametric(function, Indices({2}), Point({L0}));

Dirac mulog_eDist(L0);
mulog_eDist.setDescription(Description(1, "mulog_e"));
ComposedDistribution::DistributionCollection eColl; eColl.add(mulog_eDist); eColl.add(Dirac(0.1)); eColl.add(Dirac(0.));
ComposedDistribution eParams(eColl);

ComposedDistribution::DistributionCollection coll;
coll.add(Beta(0.117284, 0.117284, 2.9, 3.1));//R
ConditionalDistribution eDist(LogNormal(L0, 0.1, 0.), eParams);

coll.add(eDist);//e
coll.add(WeibullMin(3.16471, 9.21097, 0.0));//p
ComposedDistribution myDistribution(coll);

Point median(dim);
for(UnsignedInteger i = 0; i < dim; ++ i)
median[i] = myDistribution.getMarginal(i).computeQuantile(0.5)[0];

/* We create a 'usual' RandomVector from the Distribution */
RandomVector vect(myDistribution);

/* We create a composite random vector */
CompositeRandomVector output(parametric, vect);

/* We create an Event from this RandomVector */
ThresholdEvent myEvent(output, Less(), 0.);

fullprint << "event=" << myEvent << std::endl;

// create the algorithm
InverseFORM algo(myEvent, "mulog_e", median);

// print result
algo.run();
const InverseFORMResult result(algo.getResult());
fullprint << "Inverse FORM result=" << result << std::endl;
assert_almost_equal(result.getParameter(), Point({-4.69056}), 0.0, 0.03);

// FORM must yield the same probability on the limit state with parameter set to the optimum
eColl[0] = Dirac(result.getParameter()[0]);
eParams = ComposedDistribution(eColl);
eDist = ConditionalDistribution(LogNormal(result.getParameter()[0], 0.1, 0.0), eParams);
coll[1] = eDist;
myDistribution = ComposedDistribution(coll);
vect = RandomVector(myDistribution);
parametric.setParameter(result.getParameter());
output = CompositeRandomVector(parametric, vect);
myEvent = ThresholdEvent(output, Less(), 0.);

Cobyla cobyla;
cobyla.setMaximumIterationNumber(100);
cobyla.setMaximumAbsoluteError(1.0e-10);
cobyla.setMaximumRelativeError(1.0e-10);
cobyla.setMaximumResidualError(1.0e-10);
cobyla.setMaximumConstraintError(1.0e-10);

FORM myFORM(cobyla, myEvent, median);
myFORM.run();
const FORMResult resultFORM(myFORM.getResult());
fullprint << "result FORM="<<resultFORM << std::endl;
assert_almost_equal(resultFORM.getHasoferReliabilityIndex(), algo.getTargetBeta(), 0.0, 0.1);
}
}
catch (TestFailed & ex)
{
std::cerr << ex << std::endl;
Expand Down
4 changes: 1 addition & 3 deletions python/src/CMakeLists.txt
Original file line number Diff line number Diff line change
Expand Up @@ -219,7 +219,5 @@ endif ()

# wheel files
configure_file (METADATA.in METADATA @ONLY)
configure_file (WHEEL.in WHEEL @ONLY)
install (DIRECTORY DESTINATION ${OTROBOPT_PYTHON_MODULE_PATH}/${PACKAGE_NAME}-${PACKAGE_VERSION}.dist-info)
install (FILES ${CMAKE_CURRENT_BINARY_DIR}/METADATA ${CMAKE_CURRENT_BINARY_DIR}/WHEEL
install (FILES ${CMAKE_CURRENT_BINARY_DIR}/METADATA
DESTINATION ${OTROBOPT_PYTHON_MODULE_PATH}/${PACKAGE_NAME}-${PACKAGE_VERSION}.dist-info)
2 changes: 1 addition & 1 deletion python/src/InverseFORM_doc.i.in
Original file line number Diff line number Diff line change
Expand Up @@ -5,7 +5,7 @@ Parameters
----------
event : :class:`~openturns.RandomVector`
Event we are computing the probability of.
parameterName : string
parameterName : str
The optimized parameter
physicalStartingPoint : sequence of float
Optimization starting point
Expand Down
4 changes: 0 additions & 4 deletions python/src/WHEEL.in

This file was deleted.