MTAnisoRoughness.cpp

Go to the documentation of this file.

#include "MTAnisoRoughness.h"
#include <gsl/gsl_math.h>

namespace gplib
  {

    MTAnisoRoughness::MTAnisoRoughness() :
      conddiffweight(1.0), anisotropyweight(1.0), strikediffweight(1.0)
      {
      }

    MTAnisoRoughness::~MTAnisoRoughness()
      {
      }

    void MTAnisoRoughness::SafeParallel(const ttranscribed &member)
      {
        const unsigned int length = member.size() / 4; //we have 4 parameters in the model, so size/4 layers
        double roughness = fabs(member(2 * length)); // the loop starts at 1, so we skip this value in the 2nd summation
        const int fitexp = GetFitExponent();
        const double anglescale = 45.;
        for (unsigned int i = 1; i < length; ++i) // for all layers except the top
          {
            roughness += conddiffweight * gsl_pow_int((member(i)
                - member(i - 1)), fitexp); //minimize difference in conductivity between layers
            roughness += anisotropyweight * gsl_pow_int(fabs(member(i + 2
                * length)), fitexp); //minimize anisotropy factor in each layer
            roughness += strikediffweight * gsl_pow_int((member(i + 3 * length)
                - member(i - 1 + 3 * length)) / anglescale, fitexp); //minimize difference in angle between layers
          }
        SetRMS(pow(roughness, 1. / fitexp));
      }

    double MTAnisoRoughness::PostParallel(const ttranscribed &member)
      {
        return GetRMS();
      }

    MTAnisoRoughness::MTAnisoRoughness(const MTAnisoRoughness &Old):
      GeneralObjective(Old), conddiffweight(Old.conddiffweight),
          anisotropyweight(Old.anisotropyweight), strikediffweight(
              Old.strikediffweight)
      {

      }
    MTAnisoRoughness& MTAnisoRoughness::operator=(
        const MTAnisoRoughness& source)
      {
        if (this == &source)
          return *this;
        GeneralObjective::operator=(source);
        conddiffweight = source.conddiffweight;
        anisotropyweight = source.anisotropyweight;
        strikediffweight = source.strikediffweight;
        return *this;
      }
  }