star-mc

commit bf4b08bacd56405870fa4a158736ab4d4ab0a51a
parent 12d96a3a84f92ed2941815a2d56dfa8cad37eb19
Author: Vincent Forest <vincent.forest@meso-star.com>
Date:   Wed,  1 Apr 2015 09:38:53 +0200

First draft of the smc_solve function

Diffstat:
M
cmake/CMakeLists.txt
 | 
2
+-
A
src/smc_integrator.c
 | 
149
+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

2 files changed, 150 insertions(+), 1 deletion(-)
diff --git a/cmake/CMakeLists.txt b/cmake/CMakeLists.txt
@@ -54,7 +54,7 @@ set(VERSION_MINOR 0)
 set(VERSION_PATCH 0)
 set(VERSION ${VERSION_MAJOR}.${VERSION_MINOR}.${VERSION_PATCH})
 
-set(SMC_FILES_SRC smc_device.c smc_type.c)
+set(SMC_FILES_SRC smc_device.c smc_integrator.c smc_type.c)
 set(SMC_FILES_INC smc.h smc_device_c.h smc_type_c.h)
 
 # Prepend each file in the `SMC_FILES_<SRC|INC>' list by `SMC_SOURCE_DIR'
diff --git a/src/smc_integrator.c b/src/smc_integrator.c
@@ -0,0 +1,149 @@
+/* Copyright (C) |Meso|Star> 2015 (contact@meso-star.com)
+ *
+ * This software is a computer program whose purpose is to manage the
+ * statistical estimation of a function.
+ *
+ * This software is governed by the CeCILL license under French law and
+ * abiding by the rules of distribution of free software. You can use,
+ * modify and/or redistribute the software under the terms of the CeCILL
+ * license as circulated by CEA, CNRS and INRIA at the following URL
+ * "http://www.cecill.info".
+ *
+ * As a counterpart to the access to the source code and rights to copy,
+ * modify and redistribute granted by the license, users are provided only
+ * with a limited warranty and the software's author, the holder of the
+ * economic rights, and the successive licensors have only limited
+ * liability.
+ *
+ * In this respect, the user's attention is drawn to the risks associated
+ * with loading, using, modifying and/or developing or reproducing the
+ * software by the user in light of its specific status of free software,
+ * that may mean that it is complicated to manipulate, and that also
+ * therefore means that it is reserved for developers and experienced
+ * professionals having in-depth computer knowledge. Users are therefore
+ * encouraged to load and test the software's suitability as regards their
+ * requirements in conditions enabling the security of their systems and/or
+ * data to be ensured and, more generally, to use and operate it in the
+ * same conditions as regards security.
+ *
+ * The fact that you are presently reading this means that you have had
+ * knowledge of the CeCILL license and that you accept its terms. */
+
+#include "smc.h"
+#include "smc_device_c.h"
+#include "smc_type_c.h"
+
+#include <rsys/mem_allocator.h>
+
+#include <limits.h>
+
+struct smc_estimator {
+  void* value;
+  void* square_value;
+  unsigned long nsamples;
+
+  struct smc_device* dev;
+  ref_T ref;
+};
+
+/*******************************************************************************
+ * Helper functions
+ ******************************************************************************/
+static void
+estimator_release(ref_T* ref)
+{
+  struct smc_estimator* estimator;
+  struct smc_device* dev;
+  ASSERT(ref);
+
+  estimator = CONTAINER_OF(ref, struct smc_estimator, ref);
+  dev = estimator->dev;
+  if(estimator->value) MEM_FREE(dev->allocator, estimator->value);
+  if(estimator->square_value) MEM_FREE(dev->allocator, estimator->value);
+  MEM_FREE(dev->allocator, estimator);
+  SMC(device_ref_put(dev));
+}
+
+/*******************************************************************************
+ * Exported functions
+ ******************************************************************************/
+res_T
+smc_solve
+  (struct smc_device* dev,
+   void (*integrand)(void* value, void* ctx),
+   const struct smc_type* type,
+   void* ctx,
+   struct smc_estimator** out_estimator)
+{
+  struct smc_estimator* estimator = NULL;
+  unsigned long i;
+  void* val = NULL;
+  res_T res = RES_OK;
+
+  if(!dev || !integrand || !type || !ctx || !out_estimator || !check_type(type)) {
+    res = RES_BAD_ARG;
+    goto error;
+  }
+
+  estimator = MEM_CALLOC(dev->allocator, 1, sizeof(struct smc_estimator));
+  if(!estimator) {
+    res = RES_MEM_ERR;
+    goto error;
+  }
+  SMC(device_ref_get(dev));
+  estimator->dev = dev;
+  ref_init(&estimator->ref);
+  estimator->value = type->create(dev->allocator);
+  if(!estimator->value) {
+    res = RES_MEM_ERR;
+    goto error;
+  }
+  estimator->square_value = type->create(dev->allocator);
+  if(!estimator->value) {
+    res = RES_MEM_ERR;
+    goto error;
+  }
+  estimator->nsamples = 0;
+
+  val = type->create(dev->allocator);
+  if(!val) {
+    res = RES_MEM_ERR;
+    goto error;
+  }
+
+  FOR_EACH(i, 0, 32) {
+    if(estimator->nsamples == ULONG_MAX) break;
+    integrand(val, ctx);
+    type->add(estimator->value, estimator->value, val);
+    type->mul(val, val, val);
+    type->add(estimator->square_value, estimator->square_value, val);
+    ++estimator->nsamples;
+  }
+
+exit:
+  if(out_estimator) *out_estimator = estimator;
+  return res;
+error:
+  if(estimator) {
+    SMC(estimator_ref_put(estimator));
+    estimator = NULL;
+  }
+  goto exit;
+}
+
+res_T
+smc_estimator_ref_get(struct smc_estimator* estimator)
+{
+  if(!estimator) return RES_BAD_ARG;
+  ref_get(&estimator->ref);
+  return RES_OK;
+}
+
+res_T
+smc_estimator_ref_put(struct smc_estimator* estimator)
+{
+  if(!estimator) return RES_BAD_ARG;
+  ref_put(&estimator->ref, estimator_release);
+  return RES_OK;
+}
+