7103969144

Committed 05 Dec 2023 04:54PM UTC coverage: 54.681% (+0.03%) from 54.655%

Build # 7103969144

Build Type

Pull #85

github

Committed by

web-flow

Commit Message

Merge c4231f49e into 92fec0984

Pull Request Pull Request #85: Xcor cmp

Run Details

7095 of 12385 branches covered (0.0%)

Branch coverage included in aggregate %.

21 of 69 new or added lines in 2 files covered. (30.43%)

2 existing lines in 1 file now uncovered.

46832 of 86236 relevant lines covered (54.31%)

2102350.33 hits per line

Source File
Press 'n' to go to next uncovered line, 'b' for previous

88.58

/numcosmo/math/ncm_fit_gsl_mms.c

/***************************************************************************
 *            ncm_fit_gsl_mms.c
 *
 *  Tue Jun 19 10:55:06 2007
 *  Copyright  2007  Sandro Dias Pinto Vitenti
 *  <vitenti@uel.br>
 ****************************************************************************/
/*
 * numcosmo
 * Copyright (C) 2012 Sandro Dias Pinto Vitenti <vitenti@uel.br>
 *
 * numcosmo is free software: you can redistribute it and/or modify it
 * under the terms of the GNU General Public License as published by the
 * Free Software Foundation, either version 3 of the License, or
 * (at your option) any later version.
 *
 * numcosmo 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 General Public License for more details.
 *
 * You should have received a copy of the GNU General Public License along
 * with this program.  If not, see <http://www.gnu.org/licenses/>.
 */

/**
 * SECTION:ncm_fit_gsl_mms
 * @title: NcmFitGSLMMS
 * @short_description: Best-fit finder -- GSL non-linear minimization (simplex) algorithms.
 *
 * This object implements a best-fit finder using the GSL non-linear
 * minimization (simplex) algorithms. It is a subclass of #NcmFit.
 *
 */

#ifdef HAVE_CONFIG_H
#  include "config.h"
#endif /* HAVE_CONFIG_H */
#include "build_cfg.h"

#include "math/ncm_fit_gsl_mms.h"
#include "math/ncm_cfg.h"
#include "ncm_enum_types.h"

#ifndef NUMCOSMO_GIR_SCAN
#include <gsl/gsl_blas.h>
#include <gsl/gsl_multimin.h>
#endif /* NUMCOSMO_GIR_SCAN */

enum
{
  PROP_0,
  PROP_ALGO,
  PROP_SIZE,
};

struct _NcmFitGSLMMS
{
  /*< private >*/
  NcmFit parent_instance;
  gsl_multimin_fminimizer *mms;
  gsl_multimin_function f;
  NcmFitGSLMMSAlgos algo;
  gchar *desc;
  NcmVector *ss;
};


G_DEFINE_TYPE (NcmFitGSLMMS, ncm_fit_gsl_mms, NCM_TYPE_FIT)

static void
ncm_fit_gsl_mms_init (NcmFitGSLMMS *fit_gsl_mms)
{
  fit_gsl_mms->mms  = NULL;
  fit_gsl_mms->algo = 0;
  fit_gsl_mms->desc = NULL;
  fit_gsl_mms->ss   = NULL;
}

static void
_ncm_fit_gsl_mms_set_property (GObject *object, guint prop_id, const GValue *value, GParamSpec *pspec)
{
  NcmFitGSLMMS *fit_gsl_mms = NCM_FIT_GSL_MMS (object);

  g_return_if_fail (NCM_IS_FIT_GSL_MMS (object));

  switch (prop_id)
  {
    case PROP_ALGO:
    {
      if (fit_gsl_mms->mms == NULL)
        fit_gsl_mms->algo = g_value_get_enum (value);
      else if (fit_gsl_mms->f.n > 0)
        ncm_fit_gsl_mms_set_algo (fit_gsl_mms, g_value_get_enum (value));

      break;
    }
    default:                                                      /* LCOV_EXCL_LINE */
      G_OBJECT_WARN_INVALID_PROPERTY_ID (object, prop_id, pspec); /* LCOV_EXCL_LINE */
      break;                                                      /* LCOV_EXCL_LINE */
  }
}

static void
_ncm_fit_gsl_mms_get_property (GObject *object, guint prop_id, GValue *value, GParamSpec *pspec)
{
  NcmFitGSLMMS *fit_gsl_mms = NCM_FIT_GSL_MMS (object);

  g_return_if_fail (NCM_IS_FIT_GSL_MMS (object));

  switch (prop_id)
  {
    case PROP_ALGO:
      g_value_set_enum (value, fit_gsl_mms->algo);
      break;
    default:                                                      /* LCOV_EXCL_LINE */
      G_OBJECT_WARN_INVALID_PROPERTY_ID (object, prop_id, pspec); /* LCOV_EXCL_LINE */
      break;                                                      /* LCOV_EXCL_LINE */
  }
}

static gdouble nc_residual_multimin_f (const gsl_vector *x, gpointer p);

static void
_ncm_fit_gsl_mms_constructed (GObject *object)
{
  /* Chain up : start */
  G_OBJECT_CLASS (ncm_fit_gsl_mms_parent_class)->constructed (object);
  {
    NcmFitGSLMMS *fit_gsl_mms = NCM_FIT_GSL_MMS (object);
    NcmFit *fit               = NCM_FIT (fit_gsl_mms);
    NcmFitState *fstate       = ncm_fit_peek_state (fit);
    NcmMSet *mset             = ncm_fit_peek_mset (fit);
    const guint fparam_len    = ncm_fit_state_get_fparam_len (fstate);
    guint i;

    if (fparam_len > 0)
    {
      fit_gsl_mms->f.f      = &nc_residual_multimin_f;
      fit_gsl_mms->f.n      = fparam_len;
      fit_gsl_mms->f.params = fit_gsl_mms;

      fit_gsl_mms->ss = ncm_vector_new (fparam_len);

      for (i = 0; i < ncm_mset_fparams_len (mset); i++)
      {
        gdouble pscale = ncm_mset_fparam_get_scale (mset, i);

        ncm_vector_set (fit_gsl_mms->ss, i, pscale * 1e-3);
      }

      ncm_fit_gsl_mms_set_algo (fit_gsl_mms, fit_gsl_mms->algo);
    }
  }
}

static void
_ncm_fit_gsl_mms_dispose (GObject *object)
{
  NcmFitGSLMMS *fit_gsl_mms = NCM_FIT_GSL_MMS (object);

  ncm_vector_clear (&fit_gsl_mms->ss);

  /* Chain up : end */
  G_OBJECT_CLASS (ncm_fit_gsl_mms_parent_class)->dispose (object);
}

static void
_ncm_fit_gsl_mms_finalize (GObject *object)
{
  NcmFitGSLMMS *fit_gsl_mms = NCM_FIT_GSL_MMS (object);

  g_clear_pointer (&fit_gsl_mms->mms, gsl_multimin_fminimizer_free);

  /* Chain up : end */
  G_OBJECT_CLASS (ncm_fit_gsl_mms_parent_class)->finalize (object);
}

static NcmFit *_ncm_fit_gsl_mms_copy_new (NcmFit *fit, NcmLikelihood *lh, NcmMSet *mset, NcmFitGradType gtype);
static void _ncm_fit_gsl_mms_reset (NcmFit *fit);
static gboolean _ncm_fit_gsl_mms_run (NcmFit *fit, NcmFitRunMsgs mtype);
static const gchar *_ncm_fit_gsl_mms_get_desc (NcmFit *fit);

static void
ncm_fit_gsl_mms_class_init (NcmFitGSLMMSClass *klass)
{
  GObjectClass *object_class = G_OBJECT_CLASS (klass);
  NcmFitClass *fit_class     = NCM_FIT_CLASS (klass);

  object_class->set_property = &_ncm_fit_gsl_mms_set_property;
  object_class->get_property = &_ncm_fit_gsl_mms_get_property;
  object_class->constructed  = &_ncm_fit_gsl_mms_constructed;
  object_class->dispose      = &_ncm_fit_gsl_mms_dispose;
  object_class->finalize     = &_ncm_fit_gsl_mms_finalize;

  g_object_class_install_property (object_class,
                                   PROP_ALGO,
                                   g_param_spec_enum ("algorithm",
                                                      NULL,
                                                      "GSL multidimensional minimization algorithm [simplex]",
                                                      NCM_TYPE_FIT_GSLMMS_ALGOS, NCM_FIT_GSL_MMS_NMSIMPLEX2,
                                                      G_PARAM_READWRITE | G_PARAM_CONSTRUCT | G_PARAM_STATIC_NAME | G_PARAM_STATIC_BLURB));

  fit_class->copy_new = &_ncm_fit_gsl_mms_copy_new;
  fit_class->reset    = &_ncm_fit_gsl_mms_reset;
  fit_class->run      = &_ncm_fit_gsl_mms_run;
  fit_class->get_desc = &_ncm_fit_gsl_mms_get_desc;
}

static NcmFit *
_ncm_fit_gsl_mms_copy_new (NcmFit *fit, NcmLikelihood *lh, NcmMSet *mset, NcmFitGradType gtype)
{
  NcmFitGSLMMS *fit_gsl_mms = NCM_FIT_GSL_MMS (fit);

  return ncm_fit_gsl_mms_new (lh, mset, gtype, fit_gsl_mms->algo);
}

static void
_ncm_fit_gsl_mms_reset (NcmFit *fit)
{
  /* Chain up : start */
  NCM_FIT_CLASS (ncm_fit_gsl_mms_parent_class)->reset (fit);
  {
    NcmFitGSLMMS *fit_gsl_mms = NCM_FIT_GSL_MMS (fit);
    NcmFitState *fstate       = ncm_fit_peek_state (fit);
    const guint fparam_len    = ncm_fit_state_get_fparam_len (fstate);

    if (fit_gsl_mms->f.n != fparam_len)
    {
      g_clear_pointer (&fit_gsl_mms->mms, gsl_multimin_fminimizer_free);
      fit_gsl_mms->f.n = fparam_len;
      fit_gsl_mms->mms = NULL;

      ncm_vector_clear (&fit_gsl_mms->ss);

      if (fparam_len > 0)
      {
        NcmMSet *mset = ncm_fit_peek_mset (fit);
        guint i;

        ncm_fit_gsl_mms_set_algo (fit_gsl_mms, fit_gsl_mms->algo);

        fit_gsl_mms->ss = ncm_vector_new (fparam_len);

        for (i = 0; i < ncm_mset_fparams_len (mset); i++)
        {
          gdouble pscale = ncm_mset_fparam_get_scale (mset, i);

          ncm_vector_set (fit_gsl_mms->ss, i, pscale * 1e-3);
        }
      }
    }
  }
}

static gboolean
_ncm_fit_gsl_mms_run (NcmFit *fit, NcmFitRunMsgs mtype)
{
  NcmFitGSLMMS *fit_gsl_mms = NCM_FIT_GSL_MMS (fit);
  NcmFitState *fstate       = ncm_fit_peek_state (fit);
  NcmMSet *mset             = ncm_fit_peek_mset (fit);
  NcmVector *fparams        = ncm_fit_state_peek_fparams (fstate);
  const gdouble prec        = ncm_fit_get_params_reltol (fit);
  const guint fparam_len    = ncm_fit_state_get_fparam_len (fstate);
  gdouble last_size         = 1e300;
  gdouble last_min          = GSL_POSINF;
  gulong still_count        = 0;
  guint restart             = 10;
  const gdouble rfac        = 0.99;
  gint status;

  if (ncm_fit_equality_constraints_len (fit) || ncm_fit_inequality_constraints_len (fit))
    g_error ("_ncm_fit_gsl_mms_run: GSL algorithms do not support constraints.");

  g_assert (fparam_len != 0);

  ncm_mset_fparams_get_vector (mset, fparams);
  gsl_multimin_fminimizer_set (fit_gsl_mms->mms, &fit_gsl_mms->f, ncm_vector_gsl (fparams), ncm_vector_gsl (fit_gsl_mms->ss));

  do {
    gdouble size;

    ncm_fit_state_add_iter (fstate, 1);
    status = gsl_multimin_fminimizer_iterate (fit_gsl_mms->mms);

    if ((ncm_fit_state_get_niter (fstate) == 1) && !gsl_finite (fit_gsl_mms->mms->fval))
    {
      ncm_fit_params_set_vector (fit, fparams);

      return FALSE;
    }

    if (status)
      if (mtype > NCM_FIT_RUN_MSGS_NONE)
        ncm_fit_log_step_error (fit, gsl_strerror (status));


    size = gsl_multimin_fminimizer_size (fit_gsl_mms->mms);

    if (size == last_size)
    {
      if (++still_count == 3)
      {
        ncm_fit_log_step_error (fit, "size do not improve [prec: %8.5e size: %8.5e]", prec, size);
        status = GSL_SUCCESS;
      }
      else
      {
        status = GSL_CONTINUE;
      }
    }
    else
    {
      still_count = 0;
      last_size   = size;
      status      = gsl_multimin_test_size (size, prec);

      if (status == GSL_ETOL)
        status = GSL_CONTINUE;
    }

    if ((restart > 0) && (status == GSL_SUCCESS))
    {
      if (fit_gsl_mms->mms->fval < (last_min * rfac))
      {
        ncm_mset_fparams_get_vector (mset, fparams);
        gsl_multimin_fminimizer_set (fit_gsl_mms->mms, &fit_gsl_mms->f, ncm_vector_gsl (fparams), ncm_vector_gsl (fit_gsl_mms->ss));
        status = GSL_CONTINUE;
        restart--;
        last_min = fit_gsl_mms->mms->fval;
      }
    }

    ncm_fit_state_set_m2lnL_curval (fstate, fit_gsl_mms->mms->fval);
    ncm_fit_log_step (fit);
  } while ((status == GSL_CONTINUE) && (ncm_fit_state_get_niter (fstate) < ncm_fit_get_maxiter (fit)));

  ncm_fit_params_set_gsl_vector (fit, fit_gsl_mms->mms->x);

  ncm_mset_fparams_get_vector (mset, fparams);
  ncm_fit_state_set_m2lnL_curval (fstate, fit_gsl_mms->mms->fval);
  ncm_fit_state_set_m2lnL_prec (fstate, last_size);

  return TRUE;
}

static gdouble
nc_residual_multimin_f (const gsl_vector *x, gpointer p)
{
  NcmFit *fit   = NCM_FIT (p);
  NcmMSet *mset = ncm_fit_peek_mset (fit);
  gdouble result;

  ncm_fit_params_set_gsl_vector (fit, x);

  if (!ncm_mset_params_valid (mset))
    return GSL_POSINF;

  ncm_fit_m2lnL_val (fit, &result);

  return result;
}

static const gchar *
_ncm_fit_gsl_mms_get_desc (NcmFit *fit)
{
  NcmFitGSLMMS *fit_gsl_mms = NCM_FIT_GSL_MMS (fit);

  if (fit_gsl_mms->desc == NULL)
    fit_gsl_mms->desc = g_strdup_printf ("GSL Multidimensional Minimization:%s",
                                         fit_gsl_mms->mms != NULL ? gsl_multimin_fminimizer_name (fit_gsl_mms->mms) :
                                         "not-set");


  return fit_gsl_mms->desc;
}

/**
 * ncm_fit_gsl_mms_new:
 * @lh: a #NcmLikelihood
 * @mset: a #NcmMSet
 * @gtype: a #NcmFitGradType
 * @algo: a #NcmFitGSLMMSAlgos
 *
 * Creates a new #NcmFitGSLMMS object from the given likelihood, model set,
 * gradient type and algorithm.
 *
 *
 * Returns: (transfer full): a new #NcmFitGSLMMS object.
 */
NcmFit *
ncm_fit_gsl_mms_new (NcmLikelihood *lh, NcmMSet *mset, NcmFitGradType gtype, NcmFitGSLMMSAlgos algo)
{
  return g_object_new (NCM_TYPE_FIT_GSL_MMS,
                       "likelihood", lh,
                       "mset", mset,
                       "grad-type", gtype,
                       "algorithm", algo,
                       NULL
                      );
}

/**
 * ncm_fit_gsl_mms_new_default:
 * @lh: a #NcmLikelihood
 * @mset: a #NcmMSet
 * @gtype: a #NcmFitGradType
 *
 * Creates a new #NcmFitGSLMMS object from the given likelihood, model set and
 * gradient type. The algorithm is set to the default one (#NCM_FIT_GSL_MMS_NMSIMPLEX2).
 *
 * Returns: (transfer full): a new #NcmFitGSLMMS object.
 */
NcmFit *
ncm_fit_gsl_mms_new_default (NcmLikelihood *lh, NcmMSet *mset, NcmFitGradType gtype)
{
  return g_object_new (NCM_TYPE_FIT_GSL_MMS,
                       "likelihood", lh,
                       "mset", mset,
                       "grad-type", gtype,
                       NULL
                      );
}

/**
 * ncm_fit_gsl_mms_new_by_name:
 * @lh: a #NcmLikelihood
 * @mset: a #NcmMSet
 * @gtype: a #NcmFitGradType
 * @algo_name: a string with the name of the algorithm to be used
 *
 * Creates a new #NcmFitGSLMMS object from the given likelihood, model set,
 * gradient type and algorithm name. If @algo_name is NULL, the default
 * algorithm (#NCM_FIT_GSL_MMS_NMSIMPLEX2) is used.
 *
 * Returns: (transfer full): a new #NcmFitGSLMMS object.
 */
NcmFit *
ncm_fit_gsl_mms_new_by_name (NcmLikelihood *lh, NcmMSet *mset, NcmFitGradType gtype, gchar *algo_name)
{
  if (algo_name != NULL)
  {
    const GEnumValue *algo = ncm_cfg_get_enum_by_id_name_nick (NCM_TYPE_FIT_GSLMMS_ALGOS,
                                                               algo_name);

    if (algo == NULL)
      g_error ("ncm_fit_gsl_mms_new_by_name: algorithm %s not found.", algo_name);

    return ncm_fit_gsl_mms_new (lh, mset, gtype, algo->value);
  }
  else
  {
    return ncm_fit_gsl_mms_new_default (lh, mset, gtype);
  }
}

/**
 * ncm_fit_gsl_mms_set_algo:
 * @fit_gsl_mms: a #NcmFitGSLMMS.
 * @algo: a #gsl_mms_algorithm.
 *
 * Sets the algorithm to be used by the #NcmFitGSLMMS object.
 *
 */
void
ncm_fit_gsl_mms_set_algo (NcmFitGSLMMS *fit_gsl_mms, NcmFitGSLMMSAlgos algo)
{
  const gsl_multimin_fminimizer_type *ncm_fit_gsl_mms_algos[] = {
    gsl_multimin_fminimizer_nmsimplex2,
    gsl_multimin_fminimizer_nmsimplex,
    gsl_multimin_fminimizer_nmsimplex2rand,
  };
  NcmFit *fit         = NCM_FIT (fit_gsl_mms);
  NcmFitState *fstate = ncm_fit_peek_state (fit);

  g_assert (fit_gsl_mms->algo < NCM_FIT_GSL_MMS_NUM_ALGOS);
  g_assert (fit_gsl_mms->f.n  > 0);

  if (fit_gsl_mms->algo != algo)
  {
    fit_gsl_mms->algo = algo;

    if (fit_gsl_mms->mms != NULL)
      gsl_multimin_fminimizer_free (fit_gsl_mms->mms);

    fit_gsl_mms->mms = NULL;

    if (fit_gsl_mms->desc != NULL)
      g_free (fit_gsl_mms->desc);

    fit_gsl_mms->desc = NULL;
  }

  if (fit_gsl_mms->mms == NULL)
    fit_gsl_mms->mms = gsl_multimin_fminimizer_alloc (ncm_fit_gsl_mms_algos[fit_gsl_mms->algo],
                                                      ncm_fit_state_get_fparam_len (fstate));
}


1	/***************************************************************************
2	* ncm_fit_gsl_mms.c
3	*
4	* Tue Jun 19 10:55:06 2007
5	* Copyright 2007 Sandro Dias Pinto Vitenti
6	* <vitenti@uel.br>
7	****************************************************************************/
8	/*
9	* numcosmo
10	* Copyright (C) 2012 Sandro Dias Pinto Vitenti <vitenti@uel.br>
11	*
12	* numcosmo is free software: you can redistribute it and/or modify it
13	* under the terms of the GNU General Public License as published by the
14	* Free Software Foundation, either version 3 of the License, or
15	* (at your option) any later version.
16	*
17	* numcosmo is distributed in the hope that it will be useful, but
18	* WITHOUT ANY WARRANTY; without even the implied warranty of
19	* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
20	* See the GNU General Public License for more details.
21	*
22	* You should have received a copy of the GNU General Public License along
23	* with this program. If not, see <http://www.gnu.org/licenses/>.
24	*/
25
26	/**
27	* SECTION:ncm_fit_gsl_mms
28	* @title: NcmFitGSLMMS
29	* @short_description: Best-fit finder -- GSL non-linear minimization (simplex) algorithms.
30	*
31	* This object implements a best-fit finder using the GSL non-linear
32	* minimization (simplex) algorithms. It is a subclass of #NcmFit.
33	*
34	*/
35
36	#ifdef HAVE_CONFIG_H
37	# include "config.h"
38	#endif /* HAVE_CONFIG_H */
39	#include "build_cfg.h"
40
41	#include "math/ncm_fit_gsl_mms.h"
42	#include "math/ncm_cfg.h"
43	#include "ncm_enum_types.h"
44
45	#ifndef NUMCOSMO_GIR_SCAN
46	#include <gsl/gsl_blas.h>
47	#include <gsl/gsl_multimin.h>
48	#endif /* NUMCOSMO_GIR_SCAN */
49
50	enum
51	{
52	PROP_0,
53	PROP_ALGO,
54	PROP_SIZE,
55	};
56
57	struct _NcmFitGSLMMS
58	{
59	/< private >/
60	NcmFit parent_instance;
61	gsl_multimin_fminimizer *mms;
62	gsl_multimin_function f;
63	NcmFitGSLMMSAlgos algo;
64	gchar *desc;
65	NcmVector *ss;
66	};
67
68
69	G_DEFINE_TYPE (NcmFitGSLMMS, ncm_fit_gsl_mms, NCM_TYPE_FIT)	2,172✔
70
71	static void
72	ncm_fit_gsl_mms_init (NcmFitGSLMMS *fit_gsl_mms)	333✔
73	{
74	fit_gsl_mms->mms = NULL;	333✔
75	fit_gsl_mms->algo = 0;	333✔
76	fit_gsl_mms->desc = NULL;	333✔
77	fit_gsl_mms->ss = NULL;	333✔
78	}	333✔
79
80	static void
81	_ncm_fit_gsl_mms_set_property (GObject object, guint prop_id, const GValue value, GParamSpec *pspec)	333✔
82	{
83	NcmFitGSLMMS *fit_gsl_mms = NCM_FIT_GSL_MMS (object);	333✔
84
85	g_return_if_fail (NCM_IS_FIT_GSL_MMS (object));	333✔
86
87	switch (prop_id)	333✔
88	{
89	case PROP_ALGO:	333✔
90	{
91	if (fit_gsl_mms->mms == NULL)	333✔
92	fit_gsl_mms->algo = g_value_get_enum (value);	333✔
93	else if (fit_gsl_mms->f.n > 0)	×
94	ncm_fit_gsl_mms_set_algo (fit_gsl_mms, g_value_get_enum (value));	×
95
96	break;
97	}
98	default: /* LCOV_EXCL_LINE */
99	G_OBJECT_WARN_INVALID_PROPERTY_ID (object, prop_id, pspec); /* LCOV_EXCL_LINE */
100	break; /* LCOV_EXCL_LINE */
101	}
102	}
103
104	static void
105	_ncm_fit_gsl_mms_get_property (GObject object, guint prop_id, GValue value, GParamSpec *pspec)	188✔
106	{
107	NcmFitGSLMMS *fit_gsl_mms = NCM_FIT_GSL_MMS (object);	188✔
108
109	g_return_if_fail (NCM_IS_FIT_GSL_MMS (object));	188✔
110
111	switch (prop_id)	188✔
112	{
113	case PROP_ALGO:	188✔
114	g_value_set_enum (value, fit_gsl_mms->algo);	188✔
115	break;	188✔
116	default: /* LCOV_EXCL_LINE */
117	G_OBJECT_WARN_INVALID_PROPERTY_ID (object, prop_id, pspec); /* LCOV_EXCL_LINE */
118	break; /* LCOV_EXCL_LINE */
119	}
120	}
121
122	static gdouble nc_residual_multimin_f (const gsl_vector *x, gpointer p);
123
124	static void
125	_ncm_fit_gsl_mms_constructed (GObject *object)	333✔
126	{
127	/* Chain up : start */
128	G_OBJECT_CLASS (ncm_fit_gsl_mms_parent_class)->constructed (object);	333✔
129	{
130	NcmFitGSLMMS *fit_gsl_mms = NCM_FIT_GSL_MMS (object);	333✔
131	NcmFit *fit = NCM_FIT (fit_gsl_mms);	333✔
132	NcmFitState *fstate = ncm_fit_peek_state (fit);	333✔
133	NcmMSet *mset = ncm_fit_peek_mset (fit);	333✔
134	const guint fparam_len = ncm_fit_state_get_fparam_len (fstate);	333✔
135	guint i;	333✔
136
137	if (fparam_len > 0)	333✔
138	{
139	fit_gsl_mms->f.f = &nc_residual_multimin_f;	327✔
140	fit_gsl_mms->f.n = fparam_len;	327✔
141	fit_gsl_mms->f.params = fit_gsl_mms;	327✔
142
143	fit_gsl_mms->ss = ncm_vector_new (fparam_len);	327✔
144
145	for (i = 0; i < ncm_mset_fparams_len (mset); i++)	1,340✔
146	{
147	gdouble pscale = ncm_mset_fparam_get_scale (mset, i);	1,013✔
148
149	ncm_vector_set (fit_gsl_mms->ss, i, pscale * 1e-3);	1,013✔
150	}
151
152	ncm_fit_gsl_mms_set_algo (fit_gsl_mms, fit_gsl_mms->algo);	327✔
153	}
154	}
155	}	333✔
156
157	static void
158	_ncm_fit_gsl_mms_dispose (GObject *object)	330✔
159	{
160	NcmFitGSLMMS *fit_gsl_mms = NCM_FIT_GSL_MMS (object);	330✔
161
162	ncm_vector_clear (&fit_gsl_mms->ss);	330✔
163
164	/* Chain up : end */
165	G_OBJECT_CLASS (ncm_fit_gsl_mms_parent_class)->dispose (object);	330✔
166	}	330✔
167
168	static void
169	_ncm_fit_gsl_mms_finalize (GObject *object)	330✔
170	{
171	NcmFitGSLMMS *fit_gsl_mms = NCM_FIT_GSL_MMS (object);	330✔
172
173	g_clear_pointer (&fit_gsl_mms->mms, gsl_multimin_fminimizer_free);	330✔
174
175	/* Chain up : end */
176	G_OBJECT_CLASS (ncm_fit_gsl_mms_parent_class)->finalize (object);	330✔
177	}	330✔
178
179	static NcmFit _ncm_fit_gsl_mms_copy_new (NcmFit fit, NcmLikelihood lh, NcmMSet mset, NcmFitGradType gtype);
180	static void _ncm_fit_gsl_mms_reset (NcmFit *fit);
181	static gboolean _ncm_fit_gsl_mms_run (NcmFit *fit, NcmFitRunMsgs mtype);
182	static const gchar _ncm_fit_gsl_mms_get_desc (NcmFit fit);
183
184	static void
185	ncm_fit_gsl_mms_class_init (NcmFitGSLMMSClass *klass)	111✔
186	{
187	GObjectClass *object_class = G_OBJECT_CLASS (klass);	111✔
188	NcmFitClass *fit_class = NCM_FIT_CLASS (klass);	111✔
189
190	object_class->set_property = &_ncm_fit_gsl_mms_set_property;	111✔
191	object_class->get_property = &_ncm_fit_gsl_mms_get_property;	111✔
192	object_class->constructed = &_ncm_fit_gsl_mms_constructed;	111✔
193	object_class->dispose = &_ncm_fit_gsl_mms_dispose;	111✔
194	object_class->finalize = &_ncm_fit_gsl_mms_finalize;	111✔
195
196	g_object_class_install_property (object_class,	111✔
197	PROP_ALGO,
198	g_param_spec_enum ("algorithm",
199	NULL,
200	"GSL multidimensional minimization algorithm [simplex]",
201	NCM_TYPE_FIT_GSLMMS_ALGOS, NCM_FIT_GSL_MMS_NMSIMPLEX2,
202	G_PARAM_READWRITE \| G_PARAM_CONSTRUCT \| G_PARAM_STATIC_NAME \| G_PARAM_STATIC_BLURB));
203
204	fit_class->copy_new = &_ncm_fit_gsl_mms_copy_new;	111✔
205	fit_class->reset = &_ncm_fit_gsl_mms_reset;	111✔
206	fit_class->run = &_ncm_fit_gsl_mms_run;	111✔
207	fit_class->get_desc = &_ncm_fit_gsl_mms_get_desc;	111✔
208	}	111✔
209
210	static NcmFit *
211	_ncm_fit_gsl_mms_copy_new (NcmFit fit, NcmLikelihood lh, NcmMSet *mset, NcmFitGradType gtype)	9✔
212	{
213	NcmFitGSLMMS *fit_gsl_mms = NCM_FIT_GSL_MMS (fit);	9✔
214
215	return ncm_fit_gsl_mms_new (lh, mset, gtype, fit_gsl_mms->algo);	9✔
216	}
217
218	static void
219	_ncm_fit_gsl_mms_reset (NcmFit *fit)	926✔
220	{
221	/* Chain up : start */
222	NCM_FIT_CLASS (ncm_fit_gsl_mms_parent_class)->reset (fit);	926✔
223	{
224	NcmFitGSLMMS *fit_gsl_mms = NCM_FIT_GSL_MMS (fit);	926✔
225	NcmFitState *fstate = ncm_fit_peek_state (fit);	926✔
226	const guint fparam_len = ncm_fit_state_get_fparam_len (fstate);	926✔
227
228	if (fit_gsl_mms->f.n != fparam_len)	926✔
229	{
230	g_clear_pointer (&fit_gsl_mms->mms, gsl_multimin_fminimizer_free);	6✔
231	fit_gsl_mms->f.n = fparam_len;	6✔
232	fit_gsl_mms->mms = NULL;	6✔
233
234	ncm_vector_clear (&fit_gsl_mms->ss);	6✔
235
236	if (fparam_len > 0)	6✔
237	{
238	NcmMSet *mset = ncm_fit_peek_mset (fit);	6✔
239	guint i;	6✔
240
241	ncm_fit_gsl_mms_set_algo (fit_gsl_mms, fit_gsl_mms->algo);	6✔
242
243	fit_gsl_mms->ss = ncm_vector_new (fparam_len);	6✔
244
245	for (i = 0; i < ncm_mset_fparams_len (mset); i++)	15✔
246	{
247	gdouble pscale = ncm_mset_fparam_get_scale (mset, i);	9✔
248
249	ncm_vector_set (fit_gsl_mms->ss, i, pscale * 1e-3);	9✔
250	}
251	}
252	}
253	}
254	}	926✔
255
256	static gboolean
257	_ncm_fit_gsl_mms_run (NcmFit *fit, NcmFitRunMsgs mtype)	917✔
258	{
259	NcmFitGSLMMS *fit_gsl_mms = NCM_FIT_GSL_MMS (fit);	917✔
260	NcmFitState *fstate = ncm_fit_peek_state (fit);	917✔
261	NcmMSet *mset = ncm_fit_peek_mset (fit);	917✔
262	NcmVector *fparams = ncm_fit_state_peek_fparams (fstate);	917✔
263	const gdouble prec = ncm_fit_get_params_reltol (fit);	917✔
264	const guint fparam_len = ncm_fit_state_get_fparam_len (fstate);	917✔
265	gdouble last_size = 1e300;	917✔
266	gdouble last_min = GSL_POSINF;	917✔
267	gulong still_count = 0;	917✔
268	guint restart = 10;	917✔
269	const gdouble rfac = 0.99;	917✔
270	gint status;	917✔
271
272	if (ncm_fit_equality_constraints_len (fit) \|\| ncm_fit_inequality_constraints_len (fit))	917✔
273	g_error ("_ncm_fit_gsl_mms_run: GSL algorithms do not support constraints.");	×
274
275	g_assert (fparam_len != 0);	917✔
276
277	ncm_mset_fparams_get_vector (mset, fparams);	917✔
278	gsl_multimin_fminimizer_set (fit_gsl_mms->mms, &fit_gsl_mms->f, ncm_vector_gsl (fparams), ncm_vector_gsl (fit_gsl_mms->ss));	917✔
279
280	do {	157,691✔
281	gdouble size;	157,691✔
282
283	ncm_fit_state_add_iter (fstate, 1);	157,691✔
284	status = gsl_multimin_fminimizer_iterate (fit_gsl_mms->mms);	157,691✔
285
286	if ((ncm_fit_state_get_niter (fstate) == 1) && !gsl_finite (fit_gsl_mms->mms->fval))	157,691✔
287	{
288	ncm_fit_params_set_vector (fit, fparams);	×
289
290	return FALSE;	×
291	}
292
293	if (status)	157,691✔
294	if (mtype > NCM_FIT_RUN_MSGS_NONE)	×
295	ncm_fit_log_step_error (fit, gsl_strerror (status));	×
296
297
298	size = gsl_multimin_fminimizer_size (fit_gsl_mms->mms);	157,691✔
299
300	if (size == last_size)	157,691✔
301	{
302	if (++still_count == 3)	202✔
303	{
UNCOV 304	ncm_fit_log_step_error (fit, "size do not improve [prec: %8.5e size: %8.5e]", prec, size);	×
UNCOV 305	status = GSL_SUCCESS;	×
306	}
307	else
308	{
309	status = GSL_CONTINUE;
310	}
311	}
312	else
313	{
314	still_count = 0;	157,489✔
315	last_size = size;	157,489✔
316	status = gsl_multimin_test_size (size, prec);	157,489✔
317
318	if (status == GSL_ETOL)	157,489✔
319	status = GSL_CONTINUE;	×
320	}
321
322	if ((restart > 0) && (status == GSL_SUCCESS))	157,691✔
323	{
324	if (fit_gsl_mms->mms->fval < (last_min * rfac))	2,465✔
325	{
326	ncm_mset_fparams_get_vector (mset, fparams);	1,626✔
327	gsl_multimin_fminimizer_set (fit_gsl_mms->mms, &fit_gsl_mms->f, ncm_vector_gsl (fparams), ncm_vector_gsl (fit_gsl_mms->ss));	1,626✔
328	status = GSL_CONTINUE;	1,626✔
329	restart--;	1,626✔
330	last_min = fit_gsl_mms->mms->fval;	1,626✔
331	}
332	}
333
334	ncm_fit_state_set_m2lnL_curval (fstate, fit_gsl_mms->mms->fval);	157,691✔
335	ncm_fit_log_step (fit);	157,691✔
336	} while ((status == GSL_CONTINUE) && (ncm_fit_state_get_niter (fstate) < ncm_fit_get_maxiter (fit)));	157,691✔
337
338	ncm_fit_params_set_gsl_vector (fit, fit_gsl_mms->mms->x);	917✔
339
340	ncm_mset_fparams_get_vector (mset, fparams);	917✔
341	ncm_fit_state_set_m2lnL_curval (fstate, fit_gsl_mms->mms->fval);	917✔
342	ncm_fit_state_set_m2lnL_prec (fstate, last_size);	917✔
343
344	return TRUE;	917✔
345	}
346
347	static gdouble
348	nc_residual_multimin_f (const gsl_vector *x, gpointer p)	237,176✔
349	{
350	NcmFit *fit = NCM_FIT (p);	237,176✔
351	NcmMSet *mset = ncm_fit_peek_mset (fit);	237,176✔
352	gdouble result;	237,176✔
353
354	ncm_fit_params_set_gsl_vector (fit, x);	237,176✔
355
356	if (!ncm_mset_params_valid (mset))	237,176✔
357	return GSL_POSINF;
358
359	ncm_fit_m2lnL_val (fit, &result);	237,176✔
360
361	return result;	237,176✔
362	}
363
364	static const gchar *
365	_ncm_fit_gsl_mms_get_desc (NcmFit *fit)	18✔
366	{
367	NcmFitGSLMMS *fit_gsl_mms = NCM_FIT_GSL_MMS (fit);	18✔
368
369	if (fit_gsl_mms->desc == NULL)	18✔
370	fit_gsl_mms->desc = g_strdup_printf ("GSL Multidimensional Minimization:%s",	9✔
371	fit_gsl_mms->mms != NULL ? gsl_multimin_fminimizer_name (fit_gsl_mms->mms) :	9✔
372	"not-set");
373
374
375	return fit_gsl_mms->desc;	18✔
376	}
377
378	/**
379	* ncm_fit_gsl_mms_new:
380	* @lh: a #NcmLikelihood
381	* @mset: a #NcmMSet
382	* @gtype: a #NcmFitGradType
383	* @algo: a #NcmFitGSLMMSAlgos
384	*
385	* Creates a new #NcmFitGSLMMS object from the given likelihood, model set,
386	* gradient type and algorithm.
387	*
388	*
389	* Returns: (transfer full): a new #NcmFitGSLMMS object.
390	*/
391	NcmFit *
392	ncm_fit_gsl_mms_new (NcmLikelihood lh, NcmMSet mset, NcmFitGradType gtype, NcmFitGSLMMSAlgos algo)	145✔
393	{
394	return g_object_new (NCM_TYPE_FIT_GSL_MMS,	145✔
395	"likelihood", lh,
396	"mset", mset,
397	"grad-type", gtype,
398	"algorithm", algo,
399	NULL
400	);
401	}
402
403	/**
404	* ncm_fit_gsl_mms_new_default:
405	* @lh: a #NcmLikelihood
406	* @mset: a #NcmMSet
407	* @gtype: a #NcmFitGradType
408	*
409	* Creates a new #NcmFitGSLMMS object from the given likelihood, model set and
410	* gradient type. The algorithm is set to the default one (#NCM_FIT_GSL_MMS_NMSIMPLEX2).
411	*
412	* Returns: (transfer full): a new #NcmFitGSLMMS object.
413	*/
414	NcmFit *
415	ncm_fit_gsl_mms_new_default (NcmLikelihood lh, NcmMSet mset, NcmFitGradType gtype)	×
416	{
417	return g_object_new (NCM_TYPE_FIT_GSL_MMS,	×
418	"likelihood", lh,
419	"mset", mset,
420	"grad-type", gtype,
421	NULL
422	);
423	}
424
425	/**
426	* ncm_fit_gsl_mms_new_by_name:
427	* @lh: a #NcmLikelihood
428	* @mset: a #NcmMSet
429	* @gtype: a #NcmFitGradType
430	* @algo_name: a string with the name of the algorithm to be used
431	*
432	* Creates a new #NcmFitGSLMMS object from the given likelihood, model set,
433	* gradient type and algorithm name. If @algo_name is NULL, the default
434	* algorithm (#NCM_FIT_GSL_MMS_NMSIMPLEX2) is used.
435	*
436	* Returns: (transfer full): a new #NcmFitGSLMMS object.
437	*/
438	NcmFit *
439	ncm_fit_gsl_mms_new_by_name (NcmLikelihood lh, NcmMSet mset, NcmFitGradType gtype, gchar *algo_name)	136✔
440	{
441	if (algo_name != NULL)	136✔
442	{
443	const GEnumValue *algo = ncm_cfg_get_enum_by_id_name_nick (NCM_TYPE_FIT_GSLMMS_ALGOS,	136✔
444	algo_name);
445
446	if (algo == NULL)	136✔
447	g_error ("ncm_fit_gsl_mms_new_by_name: algorithm %s not found.", algo_name);	×
448
449	return ncm_fit_gsl_mms_new (lh, mset, gtype, algo->value);	136✔
450	}
451	else
452	{
453	return ncm_fit_gsl_mms_new_default (lh, mset, gtype);	×
454	}
455	}
456
457	/**
458	* ncm_fit_gsl_mms_set_algo:
459	* @fit_gsl_mms: a #NcmFitGSLMMS.
460	* @algo: a #gsl_mms_algorithm.
461	*
462	* Sets the algorithm to be used by the #NcmFitGSLMMS object.
463	*
464	*/
465	void
466	ncm_fit_gsl_mms_set_algo (NcmFitGSLMMS *fit_gsl_mms, NcmFitGSLMMSAlgos algo)	333✔
467	{
468	const gsl_multimin_fminimizer_type *ncm_fit_gsl_mms_algos[] = {	333✔
469	gsl_multimin_fminimizer_nmsimplex2,
470	gsl_multimin_fminimizer_nmsimplex,
471	gsl_multimin_fminimizer_nmsimplex2rand,
472	};
473	NcmFit *fit = NCM_FIT (fit_gsl_mms);	333✔
474	NcmFitState *fstate = ncm_fit_peek_state (fit);	333✔
475
476	g_assert (fit_gsl_mms->algo < NCM_FIT_GSL_MMS_NUM_ALGOS);	333✔
477	g_assert (fit_gsl_mms->f.n > 0);	333✔
478
479	if (fit_gsl_mms->algo != algo)	333✔
480	{
481	fit_gsl_mms->algo = algo;	×
482
483	if (fit_gsl_mms->mms != NULL)	×
484	gsl_multimin_fminimizer_free (fit_gsl_mms->mms);	×
485
486	fit_gsl_mms->mms = NULL;	×
487
488	if (fit_gsl_mms->desc != NULL)	×
489	g_free (fit_gsl_mms->desc);	×
490
491	fit_gsl_mms->desc = NULL;	×
492	}
493
494	if (fit_gsl_mms->mms == NULL)	333✔
495	fit_gsl_mms->mms = gsl_multimin_fminimizer_alloc (ncm_fit_gsl_mms_algos[fit_gsl_mms->algo],	333✔
496	ncm_fit_state_get_fparam_len (fstate));	333✔
497	}	333✔
498

NumCosmo / NumCosmo / 7103969144

Source File Press 'n' to go to next uncovered line, 'b' for previous

Source File
Press 'n' to go to next uncovered line, 'b' for previous