Commit 7822b856 authored by Thom Badings's avatar Thom Badings
Browse files

Push

parent a428a999
......@@ -20,9 +20,15 @@ python3 samplingprism.py --model-file crowds/crowds_10_5.pm --property-file crow
python3 samplingprism.py --model-file crowds/crowds_10_5.pm --property-file crowds/crowds.prctl --threshold 0.9 --bisimulation weak --direction geq --num_samples '(100,500,1000,5000,10000)' --num_iter 10 --beta '(0.9,0.99,0.999,0.9999)' --eta 0.413 --outfile crows_10-5_unsat_output.txt
### CROWDS (15,7)
python3 samplingprism.py --model-file crowds/crowds_15_7.pm --property-file crowds/crowds.prctl --threshold 0.9 --bisimulation strong --direction leq --num_samples '(100,500,1000,5000,10000)' --num_iter 10 --beta '(0.9,0.99,0.999,0.9999)' --eta 0.411 --outfile crowds_15-7_sat_output.txt
python3 samplingprism.py --model-file crowds/crowds_15_7.pm --property-file crowds/crowds.prctl --threshold 0.9 --bisimulation strong --direction geq --num_samples '(100,500,1000,5000,10000)' --num_iter 10 --beta '(0.9,0.99,0.999,0.9999)' --eta 0.539 --outfile crowds_15-7_unsat_output.txt
### CROWDS (20,7)
python3 samplingprism.py --model-file crowds/crowds_20_7.pm --property-file crowds/crowds.prctl --threshold 0.9 --bisimulation strong --direction leq --num_samples '(100,500,1000,5000,10000)' --num_iter 5 --beta '(0.9,0.99,0.999,0.9999)' --eta 0.416 --outfile crowds_20-7_sat_output.txt --verbose 1
python3 samplingprism.py --model-file crowds/crowds_20_7.pm --property-file crowds/crowds.prctl --threshold 0.9 --bisimulation strong --direction leq --num_samples '(100,500,1000,5000,10000)' --num_iter 10 --beta '(0.9,0.99,0.999,0.9999)' --eta 0.416 --outfile crowds_20-7_sat_output.txt --verbose 1
python3 samplingprism.py --model-file crowds/crowds_20_7.pm --property-file crowds/crowds.prctl --threshold 0.9 --bisimulation strong --direction geq --num_samples '(100,500,1000,5000,10000)' --num_iter 5 --beta '(0.9,0.99,0.999,0.9999)' --eta 0.534 --outfile crowds_20-7_unsat_output.txt
......
......@@ -165,18 +165,18 @@ def create_command(c,w,x,y,z):
commands = [create_command(c,w,x,y,z) for (c,w,x,y,z) in itertools.product(conditions, winddirs, xzones, yzones, zzones)]
def create_reward(w,x,y,z):
guard = "wdir = {} & xarea{} & yarea{} & zarea{} & valid".format(w, x, y, z)
#print(guard)
#updates = [get_update(w, i) for i in range(7)]
#probs = ["pc{}w{}x{}y{}z{}eff{}".format(c, w, x, y, z, i) for i in range(7)]
rewards = "RewardEnvw{}x{}y{}z{}".format(w, x, y, z,)
#print(rewards)
#weighted_updates = ["{}".format(p,) for p in zip(rewards)]
#update_str = " + ".join(weighted_updates)
return "\t[env] true & {}: {};".format(guard, rewards)
rewards = [create_reward(w,x,y,z) for (w,x,y,z) in itertools.product(winddirs, xzones, yzones, zzones)]
# def create_reward(w,x,y,z):
# guard = "wdir = {} & xarea{} & yarea{} & zarea{} & valid".format(w, x, y, z)
# #print(guard)
# #updates = [get_update(w, i) for i in range(7)]
# #probs = ["pc{}w{}x{}y{}z{}eff{}".format(c, w, x, y, z, i) for i in range(7)]
# rewards = "RewardEnvw{}x{}y{}z{}".format(w, x, y, z,)
# #print(rewards)
# #weighted_updates = ["{}".format(p,) for p in zip(rewards)]
# #update_str = " + ".join(weighted_updates)
# return "\t[env] true & {}: {};".format(guard, rewards)
# rewards = [create_reward(w,x,y,z) for (w,x,y,z) in itertools.product(winddirs, xzones, yzones, zzones)]
def parameter_groups():
groups = []
......@@ -187,7 +187,7 @@ def parameter_groups():
for (c,w,x,y,z) in itertools.product(conditions, winddirs, xzones, yzones, zzones):
groups.append(["pc{}w{}x{}y{}z{}eff{}".format(c, w, x, y, z, i) for i in range(7)])
#for r in range(1):
groups.append(["RewardEnv","RewardChd"])
# groups.append(["RewardEnv","RewardChd"])
#for (w) in (winddirs):
# groups.append(["RewardEnvw{}x{}y{}z{}".format(w,x,y,z) for x,y,z in itertools.product(xzones,yzones,zzones)])
#print(groups)
......@@ -245,28 +245,28 @@ def parameter_dirichlet_samples(parameter_instantiations):
return instantiations_sample
groups = parameter_groups()
# groups = parameter_groups()
parameter_defs = parameter_definitions(groups,False)
parameter_instantiations=parameter_dirichlet_instantiations(groups)
parameter_samples=parameter_dirichlet_samples(parameter_instantiations)
# parameter_defs = parameter_definitions(groups,False)
# parameter_instantiations=parameter_dirichlet_instantiations(groups)
# parameter_samples=parameter_dirichlet_samples(parameter_instantiations)
for element in parameter_samples:
print(np.sum(parameter_samples[tuple(element)]))
print(parameter_instantiations)
print(parameter_defs)
# for element in parameter_samples:
# print(np.sum(parameter_samples[tuple(element)]))
# print(parameter_instantiations)
# print(parameter_defs)
with open("drone_par_rewardparam.nm",'w') as file:
file.write(template_header)
file.write(area_definitions)
file.write("\n\n")
file.write(parameter_defs)
file.write(control_module)
file.write(mod_head)
file.write("\n".join(commands))
file.write(mod_foot)
file.write(reward_head)
file.write("\n".join(rewards))
file.write(reward_foot)
# with open("drone_par_param_example.nm",'w') as file:
# file.write(template_header)
# file.write(area_definitions)
# file.write("\n\n")
# file.write(parameter_defs)
# file.write(control_module)
# file.write(mod_head)
# file.write("\n".join(commands))
# file.write(mod_foot)
# file.write(reward_head)
# file.write("\n".join(rewards))
# file.write(reward_foot)
(stormpyenv-upt) mcubuktepe@asg-a32265:~/research/codes/scenario_code/drone$ python3 samplingprism_drone_scenario.py --model-file drone_par.nm --property-file drone_par.prctl --threshold 0.9 --bisimulation none --direction geq --num_samples 1000 --num_iter 5 --violation_prob 0.563
Building model from drone_par.nm
ModelType.MDP
Model supports parameters: True
Number of states before bisim: 23513
Number of params before bisim: 932
Number of states after bisim: 23513
Number of params after bisim: 932
40%|██████████████████ | 2/5 [02:03<03:00, 60.29s/it]Average violation so far: 0.398
60%|███████████████████████████ | 3/5 [03:04<02:01, 60.72s/it]Average violation so far: 0.39666666666666667
80%|████████████████████████████████████ | 4/5 [04:12<01:02, 62.88s/it]Average violation so far: 0.39725
100%|█████████████████████████████████████████████| 5/5 [05:20<00:00, 64.10s/it]
Average solver time : 64.10429401397705
printing specification violation array
[0.392, 0.404, 0.394, 0.399, 0.415]
probability of violating the spec less then the threshold for each iter:
[0.00227671 0.01913301 0.00336893 0.00840672 0.08523369]
average alpha_nu values of the array:
0.023683812925771354
......@@ -2775,25 +2775,25 @@ module mod
[chd] cond = 0 & wdir=0 -> c0pdircm2: (wdir' = 6) + c0pdircm1: (wdir' = 7) + c0pdirc0: (wdir' = wdir) + c0pdircp1: (wdir' = 1) + c0pdircp2: (wdir' = 2);
[chd] cond = 0 & wdir=1 -> c0pdircm2: (wdir' = 7) + c0pdircm1: (wdir' = 0) + c0pdirc0: (wdir' = wdir) + c0pdircp1: (wdir' = 2) + c0pdircp2: (wdir' = 3);
[chd] cond = 0 & wdir=2 -> c0pdircm2: (wdir' = 0) + c0pdircm1: (wdir' = 1) + c0pdirc0: (wdir' = wdir) + c0pdircp1: (wdir' = 3) + c0pdircp2: (wdir' = 4);
[chd] cond = 0 & (wdir>=3 & wdir<=5) -> c0pdircm2: (wdir' = wdir-2) + c0pdircm1: (wdir' = wdir-1) + c0pdirc0: (wdir' = wdir) + c0pdircp1: (wdir' = wdir+1) + c0pdircp2: (wdir' = wdir+2);
[chd] cond = 0 & (wdir>=3 | wdir<=5) -> c0pdircm2: (wdir' = wdir-2) + c0pdircm1: (wdir' = wdir-1) + c0pdirc0: (wdir' = wdir) + c0pdircp1: (wdir' = wdir+1) + c0pdircp2: (wdir' = wdir+2);
[chd] cond = 0 & wdir=6 -> c0pdircm2: (wdir' = 4) + c0pdircm1: (wdir' = 5) + c0pdirc0: (wdir' = wdir) + c0pdircp1: (wdir' = 7) + c0pdircp2: (wdir' = 0);
[chd] cond = 0 & wdir=7 -> c0pdircm2: (wdir' = 5) + c0pdircm1: (wdir' = 6) + c0pdirc0: (wdir' = wdir) + c0pdircp1: (wdir' = 0) + c0pdircp2: (wdir' = 1);
[chd] cond = 1 & wdir=0 -> c1pdircm2: (wdir' = 6) + c1pdircm1: (wdir' = 7) + c1pdirc0: (wdir' = wdir) + c1pdircp1: (wdir' = 1) + c1pdircp2: (wdir' = 2);
[chd] cond = 1 & wdir=1 -> c1pdircm2: (wdir' = 7) + c1pdircm1: (wdir' = 0) + c1pdirc0: (wdir' = wdir) + c1pdircp1: (wdir' = 2) + c1pdircp2: (wdir' = 3);
[chd] cond = 1 & wdir=2 -> c1pdircm2: (wdir' = 0) + c1pdircm1: (wdir' = 1) + c1pdirc0: (wdir' = wdir) + c1pdircp1: (wdir' = 3) + c1pdircp2: (wdir' = 4);
[chd] cond = 1 & (wdir>=3 & wdir<=5) -> c1pdircm2: (wdir' = wdir-2) + c1pdircm1: (wdir' = wdir-1) + c1pdirc0: (wdir' = wdir) + c1pdircp1: (wdir' = wdir+1) + c1pdircp2: (wdir' = wdir+2);
[chd] cond = 1 & (wdir>=3 | wdir<=5) -> c1pdircm2: (wdir' = wdir-2) + c1pdircm1: (wdir' = wdir-1) + c1pdirc0: (wdir' = wdir) + c1pdircp1: (wdir' = wdir+1) + c1pdircp2: (wdir' = wdir+2);
[chd] cond = 1 & wdir=6 -> c1pdircm2: (wdir' = 4) + c1pdircm1: (wdir' = 5) + c1pdirc0: (wdir' = wdir) + c1pdircp1: (wdir' = 7) + c1pdircp2: (wdir' = 0);
[chd] cond = 1 & wdir=7 -> c1pdircm2: (wdir' = 5) + c1pdircm1: (wdir' = 6) + c1pdirc0: (wdir' = wdir) + c1pdircp1: (wdir' = 0) + c1pdircp2: (wdir' = 1);
[chd] cond = 2 & wdir=0 -> c2pdircm2: (wdir' = 6) + c2pdircm1: (wdir' = 7) + c2pdirc0: (wdir' = wdir) + c2pdircp1: (wdir' = 1) + c2pdircp2: (wdir' = 2);
[chd] cond = 2 & wdir=1 -> c2pdircm2: (wdir' = 7) + c2pdircm1: (wdir' = 0) + c2pdirc0: (wdir' = wdir) + c2pdircp1: (wdir' = 2) + c2pdircp2: (wdir' = 3);
[chd] cond = 2 & wdir=2 -> c2pdircm2: (wdir' = 0) + c2pdircm1: (wdir' = 1) + c2pdirc0: (wdir' = wdir) + c2pdircp1: (wdir' = 3) + c2pdircp2: (wdir' = 4);
[chd] cond = 2 & (wdir>=3 & wdir<=5) -> c2pdircm2: (wdir' = wdir-2) + c2pdircm1: (wdir' = wdir-1) + c2pdirc0: (wdir' = wdir) + c2pdircp1: (wdir' = wdir+1) + c2pdircp2: (wdir' = wdir+2);
[chd] cond = 2 & (wdir>=3 | wdir<=5) -> c2pdircm2: (wdir' = wdir-2) + c2pdircm1: (wdir' = wdir-1) + c2pdirc0: (wdir' = wdir) + c2pdircp1: (wdir' = wdir+1) + c2pdircp2: (wdir' = wdir+2);
[chd] cond = 2 & wdir=6 -> c2pdircm2: (wdir' = 4) + c2pdircm1: (wdir' = 5) + c2pdirc0: (wdir' = wdir) + c2pdircp1: (wdir' = 7) + c2pdircp2: (wdir' = 0);
[chd] cond = 2 & wdir=7 -> c2pdircm2: (wdir' = 5) + c2pdircm1: (wdir' = 6) + c2pdirc0: (wdir' = wdir) + c2pdircp1: (wdir' = 0) + c2pdircp2: (wdir' = 1);
[chd] cond = 3 & wdir=0 -> c3pdircm2: (wdir' = 6) + c3pdircm1: (wdir' = 7) + c3pdirc0: (wdir' = wdir) + c3pdircp1: (wdir' = 1) + c3pdircp2: (wdir' = 2);
[chd] cond = 3 & wdir=1 -> c3pdircm2: (wdir' = 7) + c3pdircm1: (wdir' = 0) + c3pdirc0: (wdir' = wdir) + c3pdircp1: (wdir' = 2) + c3pdircp2: (wdir' = 3);
[chd] cond = 3 & wdir=2 -> c3pdircm2: (wdir' = 0) + c3pdircm1: (wdir' = 1) + c3pdirc0: (wdir' = wdir) + c3pdircp1: (wdir' = 3) + c3pdircp2: (wdir' = 4);
[chd] cond = 3 & (wdir>=3 & wdir<=5) -> c3pdircm2: (wdir' = wdir-2) + c3pdircm1: (wdir' = wdir-1) + c3pdirc0: (wdir' = wdir) + c3pdircp1: (wdir' = wdir+1) + c3pdircp2: (wdir' = wdir+2);
[chd] cond = 3 & (wdir>=3 | wdir<=5) -> c3pdircm2: (wdir' = wdir-2) + c3pdircm1: (wdir' = wdir-1) + c3pdirc0: (wdir' = wdir) + c3pdircp1: (wdir' = wdir+1) + c3pdircp2: (wdir' = wdir+2);
[chd] cond = 3 & wdir=6 -> c3pdircm2: (wdir' = 4) + c3pdircm1: (wdir' = 5) + c3pdirc0: (wdir' = wdir) + c3pdircp1: (wdir' = 7) + c3pdircp2: (wdir' = 0);
[chd] cond = 3 & wdir=7 -> c3pdircm2: (wdir' = 5) + c3pdircm1: (wdir' = 6) + c3pdirc0: (wdir' = wdir) + c3pdircp1: (wdir' = 0) + c3pdircp2: (wdir' = 1);
......
This diff is collapsed.
This source diff could not be displayed because it is too large. You can view the blob instead.
File mode changed from 100644 to 100755
from __future__ import division
import stormpy
import stormpy.core
import stormpy.logic
import stormpy.pars
import re
import stormpy.examples
import stormpy.examples.files
import time
import numpy as np
import time
from tqdm import tqdm
import click
import create_prism
from gurobipy import *
import argparse
# Define the parser
parser = argparse.ArgumentParser(description='Short sample app')
# Declare an argument (`--algo`), saying that the
# corresponding value should be stored in the `algo`
# field, and using a default value if the argument
# isn't given
parser.add_argument('--model-file', action="store", dest='model', default=0)
parser.add_argument('--property-file', action="store", dest='property', default=0)
parser.add_argument('--bisimulation', action="store", dest='bisim', default=0)
parser.add_argument('--threshold', action="store", dest='threshold', default=0)
parser.add_argument('--num_samples', action="store", dest='num_samples', default=0)
parser.add_argument('--num_iter', action="store", dest='num_iter', default=0)
parser.add_argument('--direction', action="store", dest='direction', default=0)
parser.add_argument('--violation_prob', action="store", dest='violation_prob', default=0)
# Now, parse the command line arguments and store the
# values in the `args` variable
args = parser.parse_args()
# Individual arguments can be accessed as attributes...
#print(args.model)
model_file=args.model
property_file=args.property
bisimulation_type=args.bisim
numiter=int(args.num_iter)
numsample=int(args.num_samples)
threshold=float(args.threshold)
violation_prob=float(args.violation_prob)
direction_type=args.direction
if direction_type=="leq":
direction=True
elif direction_type=="geq":
direction=False
else:
raise RuntimeError("Invalid direction type, direction type should be 'leq' or 'geq'")
# Approximation of binomial cdf with continuity correction for large n
# n: trials, p: success prob, m: starting successes
def BCDF(p, n, m):
return 1-CDF((m-0.5-(n*p))/math.sqrt(n*p*(1-p)))
def CDF(x):
return (1.0 + math.erf(x/math.sqrt(2.0)))/2.0
def search(values, searchFor):
for k in values:
# for v in values[k]:
if searchFor in k:
return k
return None
# paraminit = [[0 for state in range(numstate)] for state in range(numstate)]
def run_sample(numiter,numsample,thres,direction,parameters,model,properties,model_file):
'''
:param numiter: number of trials to compute the number of samples that satisfies the probability
:param numsample: number of sampled pMDPs/pMCs to check
:param thres: specification threshold
:param direction: if True, then the spec is \leq, if False, then it is \geq
:return:
'''
#storing the approximate satisfaction probability for each iteration
counterarray= [0 for _ in range(numiter)]
if str(model.model_type)=='ModelType.DTMC':
instantiator = stormpy.pars.PDtmcInstantiator(model)
elif str(model.model_type)=='ModelType.MDP':
instantiator = stormpy.pars.PMdpInstantiator(model)
else:
raise RuntimeError("Invalid model type, model should be pDTMC or pMDP.")
groups = create_prism.parameter_groups()
print(groups)
#print(groups)
#instantiate parameters
parameter_defs = create_prism.parameter_definitions(groups,False)
parameter_instantiations=create_prism.parameter_dirichlet_instantiations(groups)
start3 = time.time()
#for each iteration
for iter in tqdm(range(numiter)):
counter=0
# for each sample
for i in range(int(numsample)):
parameter_samples=create_prism.parameter_dirichlet_samples(parameter_instantiations)
# print(parameter_samples)
point = dict()
for x in parameters:
#instantiate parameters
parameter_group = search(parameter_instantiations, str(x.name))
element_int=0
for element in parameter_group:
parameter_sample_array = parameter_samples[tuple(parameter_group)]
#print(parameter_sample_array)
if (str(element)) == (str(x.name)):
point[x] = parameter_sample_array[0,element_int]
#print(point[x])
element_int=element_int+1
rational_parameter_assignments = dict(
[[x, stormpy.RationalRF(val)] for x, val in point.items()])
instantiated_model = instantiator.instantiate(rational_parameter_assignments)
result = stormpy.model_checking(instantiated_model, properties[0]).at(instantiated_model.initial_states[0])
#print(result)
#append the counter according to the spec
if direction==True:
if float(result)>thres:
counter=counter+1
else:
if float(result)<thres:
counter=counter+1
counterarray[iter]=counter/((numsample)*1.0)
if iter>1:
#print(counterarray[0:iter])
print("Average violation so far:",np.mean(counterarray[0:iter]))
#print(dict1)
t3 = time.time()
print("Average solver time :", (t3 - start3)/(numiter))
return counterarray
def load_problem(model_file, property_file,bisimulation_type):
path = model_file
prism_program = stormpy.parse_prism_program(path)
print("Building model from {}".format(path))
with open(property_file) as f:
formula_str= " ".join([l.rstrip() for l in f])
#print(formula_str)
properties = stormpy.parse_properties_for_prism_program(formula_str, prism_program)
model = stormpy.build_parametric_model(prism_program, properties)
print((model.model_type))
print("Model supports parameters: {}".format(model.supports_parameters))
parameters = model.collect_probability_parameters()
numstate=model.nr_states
print("Number of states before bisim:",numstate)
# #assert len(parameters) == 2
print ("Number of params before bisim:",len(parameters))
# print(model.model_type)
# instantiator = stormpy.pars.PDtmcInstantiator(model)
#print (model.initial_states)
#gathering parameters in the bisimulated mpdel
if bisimulation_type=='strong':
model = stormpy.perform_bisimulation(model, properties, stormpy.BisimulationType.STRONG)
elif bisimulation_type=='weak':
model = stormpy.perform_bisimulation(model, properties, stormpy.BisimulationType.WEAK)
elif bisimulation_type=='none':
pass
else:
raise RuntimeError("invalid bisimulation type, bisimulation type should be 'strong', 'weak', or 'none'. ")
parameters= model.collect_probability_parameters()
parameters_rew = model.collect_reward_parameters()
parameters.update(parameters_rew)
numstate=model.nr_states
print("Number of states after bisim:",numstate)
# #assert len(parameters) == 2
print ("Number of params after bisim:",len(parameters))
return parameters,model,properties
def compute_avg_satprob(counterarray,N,eps,flag):
'''
:param counterarray: approximate satisfaction probs
:param N: number of samples
:param eps:
:param direction: if True, then the spec is \geq, if False, then it is \leq
:return:
'''
#storing probabilities for each iteration
valuearray = np.zeros(len(counterarray))
for iters in range(len(counterarray)):
#compute number of constraints to remove in the LP
if flag:
removeconstraints = int(N * (counterarray[iters]))
else:
removeconstraints = int(N * (1 - counterarray[iters]))
# print(N,removeconstraints,eps)
#approximately compute the satisfaction prob
val3 = BCDF(eps, N, removeconstraints)
valuearray[iters] = val3
if flag:
print("probability of violating the spec less then the threshold for each iter:")
print(1-valuearray)
print("average alpha_nu values of the array:")
print(1-np.mean(valuearray))
else:
print("probability of violating the spec less then the threshold for each iter:")
print(valuearray)
print("average alpha_nu values of the array:")
print(np.mean(valuearray))
#parameters,model,properties=loader()
parameters,model,properties=load_problem(model_file,property_file,bisimulation_type)
counter_array=run_sample(numiter,numsample,threshold,direction,parameters,model,properties,model_file)
print("printing specification violation array")
print(counter_array)
#this is violation probability
compute_avg_satprob(counter_array,numsample,violation_prob,direction)
import itertools
import math
import numpy as np
template_header = """
mdp
const int MAXX = 10;
const int MAXY = 10;
const int MAXZ = 10;
formula valid = (1 < x & x < MAXX-1 & 1 < y & y < MAXY-1 & 1 < z & z < MAXZ-1) & ok=1;
formula attarget = x > MAXX - 3;
formula attarget2 = x > MAXX - 1 & y > MAXY - 2 & z > MAXZ - 3;
"""
control_module = """
module control
c : [0..3] init 0;
[move] c = 0 -> 1: (c'=1);
[env] c = 1 -> 1: (c'=2);
[chc] c = 2 -> 1: (c'=3);
[chd] c = 3 -> 1: (c'=0);
endmodule
"""
mod_head = """
module mod
ok : [0..1] init 1;
cond : [0..3] init 0;
wdir : [0..7] init 0;
x : [0..MAXX] init 2;
y : [0..MAXY] init floor(MAXY/2);
z : [0..MAXZ] init floor(MAXZ/2);
phase : [0..1] init 0;
[move] x < MAXX -> 1: (x'=x+1);
[move] x > 0 -> 1: (x'=x-1);
[move] y < MAXY -> 1: (y'=y+1);
[move] y > 0 -> 1: (y'=y-1);
[move] z < MAXX -> 1: (z'=z+1);
[move] z > 0 -> 1: (z'=z-1);
[move] attarget -> 1: (phase'=1);
[chc] cond = 0 -> p00: (cond' = 0) + p01: (cond' = 1) + p02: (cond' = 2) + p03: (cond' = 3);
[chc] cond = 1 -> p10: (cond' = 0) + p11: (cond' = 1) + p12: (cond' = 2) + p13: (cond' = 3);
[chc] cond = 2 -> p20: (cond' = 0) + p21: (cond' = 1) + p22: (cond' = 2) + p23: (cond' = 3);
[chc] cond = 3 -> p30: (cond' = 0) + p31: (cond' = 1) + p32: (cond' = 2) + p33: (cond' = 3);
[chd] cond = 0 -> c0pdircm2: (wdir' = mod(wdir + 6,8)) + c0pdircm1: (wdir' = mod(wdir + 7,8)) + c0pdirc0: (wdir' = wdir) + c0pdircp1: (wdir' = mod(wdir + 1,8)) + c0pdircp2: (wdir' = mod(wdir + 2,8));
[chd] cond = 1 -> c1pdircm2: (wdir' = mod(wdir + 6,8)) + c1pdircm1: (wdir' = mod(wdir + 7,8)) + c1pdirc0: (wdir' = wdir) + c1pdircp1: (wdir' = mod(wdir + 1,8)) + c1pdircp2: (wdir' = mod(wdir + 2,8));
[chd] cond = 2 -> c2pdircm2: (wdir' = mod(wdir + 6,8)) + c2pdircm1: (wdir' = mod(wdir + 7,8)) + c2pdirc0: (wdir' = wdir) + c2pdircp1: (wdir' = mod(wdir + 1,8)) + c2pdircp2: (wdir' = mod(wdir + 2,8));
[chd] cond = 3 -> c3pdircm2: (wdir' = mod(wdir + 6,8)) + c3pdircm1: (wdir' = mod(wdir + 7,8)) + c3pdirc0: (wdir' = wdir) + c3pdircp1: (wdir' = mod(wdir + 1,8)) + c3pdircp2: (wdir' = mod(wdir + 2,8));
[env] !valid -> 1: (ok'=0) & (x'=0) & (y'=0) & (z'=0) & (phase'=0);
"""
mod_foot="""
endmodule
"""
reward_head="""
rewards "time" // flight time
"""
reward_foot="""
endrewards
"""
labels="""
"""
#import itertools
#for element in itertools.product(*somelists):
# print(element)
xsplits = [4,6]
ysplits = [7]
zsplits = [5]
x_area_splits = [0] + xsplits + ["MAXX"]
y_area_splits = [0] + ysplits + ["MAXY"]
z_area_splits = [0] + zsplits + ["MAXZ"]
def _areadefinitions(dir, splits):
return "\n".join(["formula {}area{} = {} <= {} & {} < {};".format(dir, i, l, dir, dir, h) for i, (l, h) in
enumerate(zip(splits[:-1], splits[1:]))])
x_area_definitions = _areadefinitions("x", x_area_splits)
y_area_definitions = _areadefinitions("y", y_area_splits)
z_area_definitions = _areadefinitions("z", z_area_splits)
area_definitions = "\n".join([x_area_definitions, y_area_definitions, z_area_definitions])
conditions = [0,1,2,3]
winddirs = [0,1,2,3,4,5,6,7]
xzones = list(range(len(xsplits)+1))
yzones = list(range(len(ysplits)+1))
zzones = list(range(len(zsplits)+1))
def get_update(dir,eff):
if dir % 2 == 0:
# straight
affected_var, other_var = ("x","y") if dir in [0,4] else ("y","x")
main_var_sym = "+" if dir < 4 else "-"
sec_var_op = (["","-1","+1"][eff%3]) if dir in [0,6] else (["","+1","-1"][eff%3])
return "({}' = {}{}{})".format(affected_var, affected_var, main_var_sym, math.ceil(1)) + " & " + "({}' = {}{})".format(other_var,other_var,sec_var_op)
else:
# diagionals
assert dir in [1,3,5,7]
xsym = "+" if dir in [1,7] else "-"
ysym = "+" if dir < 4 else "-"
if eff == 0:
return "(x'=x) & (y'=y)"
if eff in [3,6]:
return "(x' = x{}{}) & (y' = y{}{})".format(xsym,1,ysym,1)
if eff == 1:
if dir == 1:
return ("(x'=x+1)")
if dir == 3:
return ("(y'=y+1)")
if dir == 5:
return ("(x'=x-1)")
if dir == 7:
return ("(y'=y-1)")
if eff == 2:
if dir == 1:
return ("(y'=y+1)")
if dir == 3:
return ("(x'=x-1)")
if dir == 5:
return ("(y'=y-1)")
if dir == 7:
return ("(x'=x+1)")
if eff == 4:
if dir == 1:
return ("(x'=x+1) & (y'=y+1)")
if dir == 3:
return ("(y'=y+1) & (x'=x-1)")
if dir == 5:
return ("(x'=x-1) & (y'=y-1)")
if dir == 7:
return ("(y'=y-1) & (x'=x+1)")
if eff == 5:
if dir == 1:
return ("(y'=y+1) & (x'=x+1)")
if dir == 3:
return ("(x'=x-1) & (y'=y+1)")
if dir == 5:
return ("(y'=y-1) & (x'=x-1)")
if dir == 7:
return ("(x'=x+1) & (y'=y-1)")
assert False
def create_command(c,w,x,y,z):
guard = "cond = {} & wdir = {} & xarea{} & yarea{} & zarea{} & valid".format(c, w, x, y, z)
updates = [get_update(w, i) for i in range(7)]
probs = ["pc{}w{}x{}y{}z{}eff{}".format(c, w, x, y, z, i) for i in range(7)]
weighted_updates = ["{}: {}".format(p, u) for p, u in zip(probs, updates)]
update_str = " + ".join(weighted_updates)
return "\t[env] {} -> {};".format(guard, update_str)
commands = [create_command(c,w,x,y,z) for (c,w,x,y,z) in itertools.product(conditions, winddirs, xzones, yzones, zzones)]
def create_reward(w,x,y,z):
guard = "wdir = {} & xarea{} & yarea{} & zarea{} & valid".format(w, x, y, z)
#print(guard)
#updates = [get_update(w, i) for i in range(7)]
#probs = ["pc{}w{}x{}y{}z{}eff{}".format(c, w, x, y, z, i) for i in range(7)]
rewards = "RewardEnvw{}x{}y{}z{}".format(w, x, y, z,)
#print(rewards)
#weighted_updates = ["{}".format(p,) for p in zip(rewards)]
#update_str = " + ".join(weighted_updates)