Skip to content
Snippets Groups Projects
Commit e44a5a08 authored by Oliver Sander's avatar Oliver Sander Committed by sander
Browse files

Remove code to measure the trust-region convergence speed

I don't think I will need that again in the foreseeable future.
If ever it should be needed again the code can be taken from the
subversion archive.

[[Imported from SVN: r10020]]
parent 93dc6f24
No related branches found
No related tags found
No related merge requests found
......@@ -44,9 +44,6 @@ mgTolerance = 1e-10
# Tolerance of the base grid solver
baseTolerance = 1e-8
# Measure convergence
instrumented = false
############################
# Problem specifications
############################
......
......@@ -279,88 +279,6 @@ int main (int argc, char *argv[]) try
vtkWriter.write(resultPath + "_" + energy + "_result");
// //////////////////////////////////////////////////////////
// Recompute and compare against exact solution
// //////////////////////////////////////////////////////////
SolutionType exactSolution = x;
// //////////////////////////////////////////////////////////////////////
// Compute mass matrix and laplace matrix to emulate L2 and H1 norms
// //////////////////////////////////////////////////////////////////////
typedef P1NodalBasis<GridType::LeafGridView,double> FEBasis;
FEBasis basis(grid->leafGridView());
OperatorAssembler<FEBasis,FEBasis> operatorAssembler(basis, basis);
LaplaceAssembler<GridType, FEBasis::LocalFiniteElement, FEBasis::LocalFiniteElement> laplaceLocalAssembler;
MassAssembler<GridType, FEBasis::LocalFiniteElement, FEBasis::LocalFiniteElement> massMatrixLocalAssembler;
typedef Dune::BCRSMatrix<Dune::FieldMatrix<double,1,1> > ScalarMatrixType;
ScalarMatrixType laplaceMatrix, massMatrix;
operatorAssembler.assemble(laplaceLocalAssembler, laplaceMatrix);
operatorAssembler.assemble(massMatrixLocalAssembler, massMatrix);
double error = std::numeric_limits<double>::max();
SolutionType intermediateSolution(x.size());
// Create statistics file
std::ofstream statisticsFile((resultPath + "trStatistics").c_str());
// Compute error of the initial iterate
typedef BlockVector<TargetSpace::CoordinateType> DifferenceType;
DifferenceType difference = computeEmbeddedDifference(exactSolution, initialIterate);
H1SemiNorm< BlockVector<TargetSpace::CoordinateType> > h1Norm(laplaceMatrix);
H1SemiNorm< BlockVector<TargetSpace::CoordinateType> > l2Norm(massMatrix);
//double oldError = std::sqrt(EnergyNorm<BCRSMatrix<FieldMatrix<double, blocksize, blocksize> >, BlockVector<FieldVector<double,blocksize> > >::normSquared(geodesicDifference, hessian));
double oldError = h1Norm(difference);
int i;
for (i=0; i<maxTrustRegionSteps; i++) {
// /////////////////////////////////////////////////////
// Read iteration from file
// /////////////////////////////////////////////////////
char iSolFilename[100];
sprintf(iSolFilename, "tmp/intermediateSolution_%04d", i);
FILE* fp = fopen(iSolFilename, "rb");
if (!fp)
DUNE_THROW(IOError, "Couldn't open intermediate solution '" << iSolFilename << "'");
for (size_t j=0; j<intermediateSolution.size(); j++) {
fread(&intermediateSolution[j], sizeof(TargetSpace), 1, fp);
}
fclose(fp);
// /////////////////////////////////////////////////////
// Compute error
// /////////////////////////////////////////////////////
difference = computeEmbeddedDifference(exactSolution, intermediateSolution);
//error = std::sqrt(EnergyNorm<BCRSMatrix<FieldMatrix<double, blocksize, blocksize> >, BlockVector<FieldVector<double,blocksize> > >::normSquared(geodesicDifference, hessian));
error = h1Norm(difference);
double convRate = error / oldError;
// Output
std::cout << "Trust-region iteration: " << i << " error : " << error << ", "
<< "convrate " << convRate << std::endl;
statisticsFile << i << " " << error << " " << convRate << std::endl;
if (error < 1e-12)
break;
oldError = error;
}
// //////////////////////////////
} catch (Exception e) {
std::cout << e << std::endl;
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment