discrete_aggregation.hpp

Go to the documentation of this file.

/*
* MAST: Multidisciplinary-design Adaptation and Sensitivity Toolkit
* Copyright (C) 2013-2020  Manav Bhatia and MAST authors
*
* This library is free software; you can redistribute it and/or
* modify it under the terms of the GNU Lesser General Public
* License as published by the Free Software Foundation; either
* version 2.1 of the License, or (at your option) any later version.
*
* This library is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
* Lesser General Public License for more details.
*
* You should have received a copy of the GNU Lesser General Public
* License along with this library; if not, write to the Free Software
* Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA
*/

#ifndef __mast_optimization_discrete_aggregation_h__
#define __mast_optimization_discrete_aggregation_h__

// MAST includes
#include <mast/base/mast_data_types.h>
#include <mast/numerics/utility.hpp>

// libMesh includes
#include <libmesh/parallel.h>

namespace MAST {
namespace Optimization {
namespace Aggregation {

template <typename ScalarType>
ScalarType
aggregate_minimum(const libMesh::Parallel::Communicator *comm,
                  const std::vector<ScalarType>         &vec,
                  const real_t                           p) {
    
    ScalarType
    v      = 0.,
    v_min  = 0.;
    
    v_min = MAST::Numerics::Utility::real_minimum(vec);
    
    if (comm) MAST::Numerics::Utility::comm_min(*comm, v_min);
    
    for (uint_t i=0; i<vec.size(); i++) {
        
        v += exp(-p * (vec[i] - v_min));
    }
    
    if (comm) MAST::Numerics::Utility::comm_sum(*comm, v);

    v = v_min - log(v) / p;
    
    return v;
}


template <typename ScalarType>
ScalarType
aggregate_minimum_sensitivity(const libMesh::Parallel::Communicator *comm,
                              const std::vector<ScalarType>         &vec,
                              const uint_t                           i,
                              const real_t                           p) {
    
    ScalarType
    v      = 0.,
    v_min  = 0.;
    
    v_min = MAST::Numerics::Utility::real_minimum(vec);

    if (comm) MAST::Numerics::Utility::comm_min(*comm, v_min);

    for (uint_t i=0; i<vec.size(); i++) {
        
        v  += exp(-p * (vec[i] - v_min));
    }
    
    if (comm) MAST::Numerics::Utility::comm_sum(*comm, v);

    v = exp(-p * (vec[i] - v_min)) / v;
    
    return v;
}



template <typename ScalarType>
ScalarType
aggregate_minimum_sensitivity(const libMesh::Parallel::Communicator *comm,
                              const std::vector<ScalarType>         &vec,
                              const std::vector<ScalarType>         &dvec,
                              const real_t                           p) {
    
    ScalarType
    dv     = 0.,
    v      = 0.,
    v_min  = 0.;
    
    v_min = MAST::Numerics::Utility::real_minimum(vec);

    if (comm) MAST::Numerics::Utility::comm_min(*comm, v_min);

    for (uint_t i=0; i<vec.size(); i++) {
        
        dv += exp(-p * (vec[i] - v_min)) * dvec[i];
        v  += exp(-p * (vec[i] - v_min));
    }
    
    if (comm) MAST::Numerics::Utility::comm_sum(*comm, dv);
    if (comm) MAST::Numerics::Utility::comm_sum(*comm, v);

    v = dv / v;
    
    return v;
}



template <typename ScalarType>
void
aggregate_minimum_denominator(const libMesh::Parallel::Communicator *comm,
                              const std::vector<ScalarType>         &vec,
                              const real_t                           p,
                              ScalarType                            &denom,
                              ScalarType                            &v_min) {
    
    denom  = 0.;
    
    v_min = MAST::Numerics::Utility::real_minimum(vec);

    if (comm) MAST::Numerics::Utility::comm_min(*comm, v_min);

    for (uint_t i=0; i<vec.size(); i++) {
        
        denom  += exp(-p * (vec[i] - v_min));
    }
    
    if (comm) MAST::Numerics::Utility::comm_sum(*comm, denom);
}



template <typename ScalarType>
ScalarType
aggregate_minimum_sensitivity(const std::vector<ScalarType>   &vec,
                              const uint_t                     i,
                              const real_t                     p,
                              const ScalarType                &denom,
                              const ScalarType                &v_min) {
    
    return exp(-p * (vec[i] - v_min)) / denom;
}


template <typename ScalarType>
ScalarType
aggregate_minimum_sensitivity(const libMesh::Parallel::Communicator *comm,
                              const std::vector<ScalarType>         &vec,
                              const std::vector<ScalarType>         &dvec,
                              const real_t                           p,
                              const ScalarType                      &denom,
                              const ScalarType                      &v_min) {
    
    ScalarType
    dv     = 0.;
    
    for (uint_t i=0; i<vec.size(); i++) {
        
        dv += exp(-p * (vec[i] - v_min)) * dvec[i];
    }
    
    if (comm) MAST::Numerics::Utility::comm_sum(*comm, dv);

   return dv / denom;
}


template <typename ScalarType>
ScalarType
aggregate_maximum(const libMesh::Parallel::Communicator *comm,
                  const std::vector<ScalarType>         &vec,
                  const real_t                           p) {
    
    ScalarType
    v      = 0.,
    v_max  = 0.;

    v_max = MAST::Numerics::Utility::real_maximum(vec);

    if (comm) MAST::Numerics::Utility::comm_max(*comm, v_max);

    for (uint_t i=0; i<vec.size(); i++) {
        
        v += exp(p * (vec[i] - v_max));
    }
    
    if (comm) MAST::Numerics::Utility::comm_sum(*comm, v);

    v = v_max + log(v) / p;
    
    return v;
}


template <typename ScalarType>
ScalarType
aggregate_maximum_sensitivity(const libMesh::Parallel::Communicator *comm,
                              const std::vector<ScalarType>         &vec,
                              const uint_t                           i,
                              const real_t                           p) {
    
    ScalarType
    v      = 0.,
    v_max  = 0.;
    
    v_max = MAST::Numerics::Utility::real_maximum(vec);

    if (comm) MAST::Numerics::Utility::comm_max(*comm, v_max);
    
    for (uint_t i=0; i<vec.size(); i++) {
        
        v += exp(p * (vec[i] - v_max));
    }
    
    if (comm) MAST::Numerics::Utility::comm_sum(*comm, v);

    v = exp(p * (vec[i] - v_max)) / v;
    
    return v;
}



template <typename ScalarType>
ScalarType
aggregate_maximum_sensitivity(const libMesh::Parallel::Communicator *comm,
                              const std::vector<ScalarType>         &vec,
                              const std::vector<ScalarType>         &dvec,
                              const real_t                           p) {
    
    ScalarType
    dv     = 0.,
    v      = 0.,
    v_max  = 0.;
    
    v_max = MAST::Numerics::Utility::real_maximum(vec);

    if (comm) MAST::Numerics::Utility::comm_max(*comm, v_max);

    for (uint_t i=0; i<vec.size(); i++) {
        
        dv += exp(p * (vec[i] - v_max)) * dvec[i];
        v  += exp(p * (vec[i] - v_max));
    }
    
    if (comm) MAST::Numerics::Utility::comm_sum(*comm, dv);
    if (comm) MAST::Numerics::Utility::comm_sum(*comm, v);

    v = dv / v;
    
    return v;
}


template <typename ScalarType>
void
aggregate_maximum_denominator(const libMesh::Parallel::Communicator *comm,
                              const std::vector<ScalarType>         &vec,
                              const real_t                           p,
                              ScalarType                            &denom,
                              ScalarType                            &v_max) {
    
    denom  = 0.;
    
    v_max = MAST::Numerics::Utility::real_maximum(vec);

    if (comm) MAST::Numerics::Utility::comm_max(*comm, v_max);
    
    for (uint_t i=0; i<vec.size(); i++) {
        
        denom  += exp(p * (vec[i] - v_max));
    }

    if (comm) MAST::Numerics::Utility::comm_sum(*comm, denom);
}



template <typename ScalarType>
ScalarType
aggregate_maximum_sensitivity(const std::vector<ScalarType>   &vec,
                              const uint_t                     i,
                              const real_t                     p,
                              const ScalarType                &denom,
                              const ScalarType                &v_max) {
    
    return exp(p * (vec[i] - v_max)) / denom;
}


template <typename ScalarType>
ScalarType
aggregate_maximum_sensitivity(const libMesh::Parallel::Communicator *comm,
                              const std::vector<ScalarType>         &vec,
                              const std::vector<ScalarType>         &dvec,
                              const real_t                           p,
                              const ScalarType                      &denom,
                              const ScalarType                      &v_max) {
    
    ScalarType
    dv     = 0.;
    
    for (uint_t i=0; i<vec.size(); i++) {
        
        dv += exp(p * (vec[i] - v_max)) * dvec[i];
    }

    if (comm) MAST::Numerics::Utility::comm_sum(*comm, dv);

    return dv / denom;
}


} // Aggregation
} // Optimization
} // MAST

#endif // __mast_optimization_discrete_aggregation_h__