Commit c3db4cd0 authored by Lambert Patrick's avatar Lambert Patrick
Browse files

with instance BinomialForSuperIndividualGen

parent da64a0f5
/**
* patrick
* @author Patrick Lambert
* @copyright Copyright (c) 2018, Irstea
*
* This program 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.
*
* This program 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/>.
*/
package miscellaneous;
import umontreal.iro.lecuyer.randvar.NormalACRGen;
import umontreal.iro.lecuyer.randvar.NormalGen;
import umontreal.iro.lecuyer.randvar.UniformGen;
import umontreal.iro.lecuyer.rng.RandomStream;
/**
*
*/
public class BinomialForSuperIndividualGen {
private long threshold;
private NormalGen normalGen;
private UniformGen uniformGen;
public BinomialForSuperIndividualGen(RandomStream randomStream, long threshold){
normalGen = new NormalACRGen(randomStream, 0., 1.); // NormalACRGen to speed up exceution
uniformGen = new UniformGen(randomStream, 0., 1.);
this.threshold = threshold;
}
public BinomialForSuperIndividualGen(RandomStream randomStream){
this(randomStream, 50L);
}
public long getSuccessNumber(long draws, double succesProbability) {
double successes = -1. ;
if (draws >= threshold) {
// approximate the binomial by a normal distribution of mean n p and variance n p (1-p)
while (successes < 0 | successes > draws) {
successes = Math.round(normalGen.nextDouble() * Math.sqrt(succesProbability * (1 - succesProbability) * draws)
+ succesProbability * draws);
}
} else {
successes = 0.;
for (long i = 0; i < draws; i++) {
if (uniformGen.nextDouble() < succesProbability) {
successes++;
}
}
}
return (long) successes;
}
}
......@@ -15,6 +15,7 @@ import environment.SeaBasin;
import environment.Time;
import environment.Time.Season;
import fr.cemagref.simaqualife.kernel.processes.AquaNismsGroupProcess;
import fr.cemagref.simaqualife.pilot.Pilot;
import miscellaneous.Miscellaneous;
import org.openide.util.lookup.ServiceProvider;
......@@ -22,6 +23,7 @@ import org.openide.util.lookup.ServiceProvider;
import com.thoughtworks.xstream.XStream;
import com.thoughtworks.xstream.io.xml.DomDriver;
import miscellaneous.BinomialForSuperIndividualGen;
import miscellaneous.Duo;
/**
......@@ -91,11 +93,25 @@ public class DisperseAndMigrateToRiverWithMultiNomDistriAndDeathBasin extends Di
*/
private long yearOfTheKillings;
/**
* the random numbers generator for binomial draws
* @unit --
*/
private transient BinomialForSuperIndividualGen aleaGen;
public static void main(String[] args) {
System.out.println((new XStream(new DomDriver()))
.toXML(new DisperseAndMigrateToRiverWithMultiNomDistriAndDeathBasin()));
}
@Override
public void initTransientParameters(Pilot pilot) {
super.initTransientParameters(pilot);
aleaGen = new BinomialForSuperIndividualGen(pilot.getRandomStream());
}
@Override
public void doProcess(DiadromousFishGroup group) {
Time time = group.getEnvironment().getTime();
......@@ -126,7 +142,8 @@ public class DisperseAndMigrateToRiverWithMultiNomDistriAndDeathBasin extends Di
if (fish.isMature()) {
// fish with homing
homingAmount = Miscellaneous.binomialForSuperIndividual(group.getPilot(), fish.getAmount(), pHoming); // seuil par d�faut fix� � 50
//homingAmount = Miscellaneous.binomialForSuperIndividual(group.getPilot(), fish.getAmount(), pHoming); // seuil par d�faut fix� � 50
homingAmount = aleaGen.getSuccessNumber(fish.getAmount(), pHoming);
// strayed fish
if (killStrayers == true && time.getYear(group.getPilot()) >= yearOfTheKillings) {
strayedAmount = 0;
......@@ -168,8 +185,8 @@ public class DisperseAndMigrateToRiverWithMultiNomDistriAndDeathBasin extends Di
RiverBasin strayerDestination = entry.getKey();
Double weight = entry.getValue();
probToGo = weight / totalWeight;
amountToGo = Miscellaneous.binomialForSuperIndividual(group.getPilot(), strayedAmount, probToGo);
//amountToGo = Miscellaneous.binomialForSuperIndividual(group.getPilot(), strayedAmount, probToGo);
amountToGo = aleaGen.getSuccessNumber(strayedAmount, probToGo);
if (amountToGo > 0) {
strayerDestination.addFish(fish.duplicateWithNewPositionAndAmount(group.getPilot(), strayerDestination, amountToGo), group);
}
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment