Laboratoire RDP » Biophysics Team » Vascular Networks in Growing Leaves

#include "compute_energy_and_forces.h"

#include "minimize_BFGS.h"

#include <nlopt.h>

#include <stdio.h>

#include <stdlib.h>

#include <cmath>

/* this function minimizes the mechanical energy */

void minimize_BFGS(int vertex_number, double *vertices, int edge_number, int *edges_vertices,

                   double *lengths, double *rest_lengths, double *edges_width, int *periodizer, double *width, double *height){

nlopt_opt opt;

opt=nlopt_create(NLOPT_LD_LBFGS,2*vertex_number+2);

nlopt_set_stopval(opt, -HUGE_VAL);

nlopt_set_maxeval(opt, -1.);

nlopt_set_maxtime(opt, -1.);

nlopt_set_xtol_rel(opt, -1.);

nlopt_set_xtol_abs1(opt, 1.e-6);

nlopt_set_ftol_rel(opt, -1.);

nlopt_set_ftol_abs(opt, -1.);

struct parameters_list {int vertex_number; int edge_number; int *edges_vertices;

                        double *lengths; double *rest_lengths; double *edges_width; int *periodizer; int anisotropy; double ref_width; double ref_height;};

double min_energy;

struct parameters_list parameters;

parameters.vertex_number=vertex_number;

parameters.edge_number=edge_number;

parameters.edges_vertices=edges_vertices;

parameters.lengths=lengths;

parameters.rest_lengths=rest_lengths;

parameters.edges_width=edges_width;

parameters.periodizer=periodizer;

parameters.anisotropy=0;

parameters.ref_width=0;

parameters.ref_height=0;

void* pparameters=&parameters;

double vertices_BFGS[2*vertex_number+2];

for (int i=0; i<2*vertex_number;i++){ vertices_BFGS[i]=vertices[i];}

vertices_BFGS[2*vertex_number]=width[0]; vertices_BFGS[2*vertex_number+1]=height[0];

/// first resolution: find the reference width and height without the anisotropic stress

int reso=nlopt_set_min_objective(opt, compute_energy_and_forces, pparameters);

nlopt_optimize(opt, vertices_BFGS, &min_energy);

parameters.anisotropy=1;

parameters.ref_width=vertices_BFGS[2*vertex_number];

parameters.ref_height=vertices_BFGS[2*vertex_number+1];

/// second resolution with the anisotropic stress

reso=nlopt_set_min_objective(opt, compute_energy_and_forces, pparameters);

nlopt_optimize(opt, vertices_BFGS, &min_energy);

for (int i=0; i<2*vertex_number;i++){vertices[i]=vertices_BFGS[i];}

width[0]=vertices_BFGS[2*vertex_number]; height[0]=vertices_BFGS[2*vertex_number+1];

nlopt_destroy(opt);

}