CARLsim  6.1.0
CARLsim: a GPU-accelerated SNN simulator
poisson_rate.cpp
Go to the documentation of this file.
1 /*
2  * Copyright (c) 2014 Regents of the University of California. All rights reserved.
3  *
4  * Redistribution and use in source and binary forms, with or without
5  * modification, are permitted provided that the following conditions
6  * are met:
7  *
8  * 1. Redistributions of source code must retain the above copyright
9  * notice, this list of conditions and the following disclaimer.
10  *
11  * 2. Redistributions in binary form must reproduce the above copyright
12  * notice, this list of conditions and the following disclaimer in the
13  * documentation and/or other materials provided with the distribution.
14  *
15  * 3. The names of its contributors may not be used to endorse or promote
16  * products derived from this software without specific prior written
17  * permission.
18  *
19  * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
20  * "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
21  * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
22  * A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR
23  * CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
24  * EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
25  * PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
26  * PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
27  * LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
28  * NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
29  * SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
30  *
31  * *********************************************************************************************** *
32  * CARLsim
33  * created by: (MDR) Micah Richert, (JN) Jayram M. Nageswaran
34  * maintained by: (MA) Mike Avery <averym@uci.edu>, (MB) Michael Beyeler <mbeyeler@uci.edu>,
35  * (KDC) Kristofor Carlson <kdcarlso@uci.edu>
36  * (TSC) Ting-Shuo Chou <tingshuc@uci.edu>
37  *
38  * CARLsim6: LN, JX, KC, KW
39  *
40  * CARLsim available from http://socsci.uci.edu/~jkrichma/CARLsim/
41  * Ver 6/13/2016
42  */
43 #include <cassert> // assert
44 #include <cstdio> // printf
45 #include <cstring> // string, memset
46 #include <cstdlib> // malloc, free, rand
47 
48 #include <poisson_rate.h>
49 #include <carlsim_definitions.h> // ALL
50 #ifndef __NO_CUDA__
51  #include <cuda_version_control.h>
52 #endif
53 
54 
55 
57 public:
58  // +++++ PUBLIC METHODS +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ //
59  Impl(int nNeur, bool onGPU): nNeur_(nNeur), onGPU_(onGPU) {
60  assert(nNeur>0);
61 
62  h_rates_ = NULL;
63  d_rates_ = NULL;
64 
65  if (onGPU) {
66 #ifndef __NO_CUDA__
67  // allocate rates on device and set to zero
68  CUDA_CHECK_ERRORS(cudaMalloc((void**)&d_rates_, sizeof(float)*nNeur));
69  CUDA_CHECK_ERRORS(cudaMemset(d_rates_, 0, sizeof(float)*nNeur));
70 #else
71  printf("Cannot use onGPU when compiled without CUDA library.\n");
72  assert(false);
73 #endif
74  } else {
75  // allocate rates on host and set to zero
76  h_rates_ = new float[nNeur];
77  memset(h_rates_, 0, sizeof(float)*nNeur);
78  }
79  }
80 
81  ~Impl() {
82  if (isOnGPU()) {
83 #ifndef __NO_CUDA__
84  // clean up device
85  if (d_rates_!=NULL) {
86  CUDA_CHECK_ERRORS(cudaFree(d_rates_)); // free memory
87  d_rates_ = NULL;
88  }
89 #endif
90  } else {
91  // clean up host
92  if (h_rates_!=NULL) {
93  delete[] h_rates_;
94  }
95  h_rates_ = NULL;
96  }
97  }
98 
99  int getNumNeurons() {
100  return nNeur_;
101  }
102 
103  // get rate of certain neuron
104  float getRate(int neurId) {
105  assert(neurId>=0 && neurId<getNumNeurons());
106 
107  if (isOnGPU()) {
108  // get data from device (might have kernel launch overhead because float is small)
109  float h_d_rate = 0.0f;
110 //Fix LN20210211
111 // warning C4715: 'PoissonRate::Impl::getRate': not all control paths return a value
112 #ifndef __NO_CUDA__
113  CUDA_CHECK_ERRORS( cudaMemcpy(&h_d_rate, &(d_rates_[neurId]), sizeof(float), cudaMemcpyDeviceToHost) );
114 #endif
115  return h_d_rate;
116  } else {
117  // data is on host
118  return h_rates_[neurId];
119  }
120  }
121 
122  // get all rates as vector
123  std::vector<float> getRates() {
124  if (isOnGPU()) {
125 //Fix LN20210211
126 // warning C4715: 'PoissonRate::Impl::getRate': not all control paths return a value
127 #ifndef __NO_CUDA__
128  // get data from device
129  float *h_d_rates = (float*)malloc(sizeof(float)*getNumNeurons());
130  CUDA_CHECK_ERRORS( cudaMemcpy(h_d_rates, d_rates_, sizeof(float)*getNumNeurons(), cudaMemcpyDeviceToHost) );
131 
132  // copy data to vec
133  std::vector<float> rates(h_d_rates, h_d_rates + getNumNeurons());
134  free(h_d_rates);
135 
136  return rates;
137 #else
138  return std::vector<float>(); // fix compiler warning C4715, this code shall be never reached
139 #endif
140  } else {
141  // data is on host
142  std::vector<float> rates(h_rates_, h_rates_ + getNumNeurons()); // copy to vec
143  return rates;
144  }
145  }
146 
147  // get pointer to rate array on CPU
148  float* getRatePtrCPU() {
149  assert(!isOnGPU());
150  return h_rates_;
151  }
152 
153  // get pointer to rate array on GPU
154  float* getRatePtrGPU() {
155  assert(isOnGPU());
156  return d_rates_;
157  }
158 
159  bool isOnGPU() {
160  return onGPU_;
161  }
162 
163  // set rate of a specific neuron
164  void setRate(int neurId, float rate) {
165  if (neurId==ALL) {
166  setRates(rate);
167  } else {
168  assert(neurId>=0 && neurId<getNumNeurons());
169  if (isOnGPU()) {
170 #ifndef __NO_CUDA__
171  // copy float to device (might have kernel launch overhead because float is small)
172  CUDA_CHECK_ERRORS( cudaMemcpy(&(d_rates_[neurId]), &rate, sizeof(float), cudaMemcpyHostToDevice) );
173 #endif
174  } else {
175  // set float in host array
176  h_rates_[neurId] = rate;
177  }
178  }
179  }
180 
181  // set rate of all neurons to same value
182  void setRates(float rate) {
183  assert(rate>=0.0f);
184 
185  std::vector<float> rates(getNumNeurons(), rate);
186  setRates(rates);
187  }
188 
189  // set rates with vector
190  void setRates(const std::vector<float>& rate) {
191  assert(rate.size()==getNumNeurons());
192 
193  if (isOnGPU()) {
194 #ifndef __NO_CUDA__
195  // copy to device
196  float *h_rates_arr = new float[getNumNeurons()];
197  std::copy(rate.begin(), rate.end(), h_rates_arr);
198  CUDA_CHECK_ERRORS( cudaMemcpy(d_rates_, h_rates_arr, sizeof(float)*getNumNeurons(), cudaMemcpyHostToDevice) );
199  delete[] h_rates_arr;
200 #endif
201  } else {
202  // set host array
203  std::copy(rate.begin(), rate.end(), h_rates_);
204  }
205  }
206 
207 
208 private:
209  // +++++ PRIVATE METHODS ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ //
210 
211  // +++++ PRIVATE STATIC PROPERTIES ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ //
212 
213  // +++++ PRIVATE PROPERTIES +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ //
214 
215  float *h_rates_;
216  float *d_rates_;
217  const int nNeur_;
218  const bool onGPU_;
219 };
220 
221 
222 // ****************************************************************************************************************** //
223 // POISSONRATE API IMPLEMENTATION
224 // ****************************************************************************************************************** //
225 
226 // create and destroy a pImpl instance
227 PoissonRate::PoissonRate(int nNeur, bool onGPU) : _impl( new Impl(nNeur, onGPU) ) {}
228 PoissonRate::~PoissonRate() { delete _impl; }
229 
230 int PoissonRate::getNumNeurons() { return _impl->getNumNeurons(); }
231 float PoissonRate::getRate(int neurId) { return _impl->getRate(neurId); }
232 std::vector<float> PoissonRate::getRates() { return _impl->getRates(); }
233 float* PoissonRate::getRatePtrCPU() { return _impl->getRatePtrCPU(); }
234 float* PoissonRate::getRatePtrGPU() { return _impl->getRatePtrGPU(); }
235 bool PoissonRate::isOnGPU() { return _impl->isOnGPU(); }
236 void PoissonRate::setRate(int neurId, float rate) { _impl->setRate(neurId, rate); }
237 void PoissonRate::setRates(float rate) { _impl->setRates(rate); }
238 void PoissonRate::setRates(const std::vector<float>& rates) { _impl->setRates(rates); }
void setRates(float rate)
std::vector< float > getRates()
void setRates(float rate)
Assigns the same mean firing rate to all neurons.
#define ALL
CARLsim common definitions.
float * getRatePtrGPU()
float * getRatePtrCPU()
Returns pointer to CPU-allocated firing rate array (deprecated)
Impl(int nNeur, bool onGPU)
float getRate(int neurId)
Returns the mean firing rate of a specific neuron ID.
float * getRatePtrCPU()
~PoissonRate()
PoissonRate destructor.
int getNumNeurons()
Returns the number of neurons for which to generate Poisson spike trains.
void setRate(int neurId, float rate)
std::vector< float > getRates()
Returns a vector of firing rates, one element per neuron.
void setRates(const std::vector< float > &rate)
void setRate(int neurId, float rate)
Sets the mean firing rate of a particular neuron ID.
PoissonRate(int nNeur, bool onGPU=false)
PoissonRate constructor.
float getRate(int neurId)
float * getRatePtrGPU()
Returns pointer to GPU-allocated firing rate array (deprecated)
bool isOnGPU()
Checks whether the firing rates are allocated on CPU or GPU.