Below is the explaination regarding this code from Jordan. Attached is the code in its last form.
Hi Lucas,
I'll try to fit a bunch of stuff into this email, so hopefully some of it is useful/interesting.
As far as the evolutionary algorithm goes, I'll send you the code and let you look at it yourself. the paper that I am developing it from comes from https://www.researchgate.net/publication/282857432_INTEGRATION_OF_REMCOM_XFDTD_AND_OCTAVE_FOR_THE_AUTOMATED_OPTIMAL_DESIGN_OF_MICROWAVE_ANTENNAS. My code has a couple differences in that I didn't look at the past k iterations like the paper does, but instead I have 5 parallel designs being run and look at how many of those improve the output. This is subtly different and might not be good, because it considers each design at the same time, so the mean doesn't have time to adjust before it is reevaluated if that makes sense. Anyways, just something to think about.
The basic structure for the code is that it is run with command line arguments, so that you compile it with whatever name ( I usually call it evolved_dipole, but doesn't matter). So it is run as
$./evolved_dipole --start
for the first time in a run and then every subsequent time in a run
$./evolve_dipole --cont
The --start will create the initial population, and record the parameters in a .csv file called handshake.csv. The --cont will read in a file called handshook.csv that theoretically will have the output of the antenna simulations for each antenna.
The first obvious thing I can think of that is missing in this script is that it doesn't write a number to a txt file in the watch folder, but I'll explain that later. The second obvious thing that I didn't add to this is the checking of the exit condition d/do < some value (see paper if this is confusing). The third thing I can think of is that I don't have any constraints on the values that the script produces. It will probably be valuable to include some constraints on these values at the very least so that you don't get negative distances. In addition, this script should be easily generalizable by increasing NVAR and changing the mean and deviation vectors. The code is also not particularly pretty so I apologize for that. I tried to comment a decent amount.
Then, there is the XF script. So the XF script should input the antenna parameters, create the simulations and then start them. Then output the data. One of the things that I never ended up doing here is that the data output can only be done after the simulations have finished running, so you'll need to figure out how to get that to work if you use XF, so I'll include the scripts separately. For the output.xmacro script you will need to specify the simulation and run id each time it is used. I think it might be possible to just have a while function that will wait x seconds and then reevaluate and the condition is whether the data is available, but that might not be the best way.
Then, we get to the controlling portion of the code. I have a batch script which should control a bash script (batch submission for clusters is weird). so the bash script theoretically should load the xf module, run evolved_dipole --start, watch a folder and see if a certain number is written to a file, run the xf simulation, watch again for a different number, run the evolved_dipole --cont (or arasim first eventually) and then it creates a loop and should run theoretically until the exit condition is reached in which case --cont will need to write a different value than before and everything should end (i don't know that I've included this part yet).
The big problem here is that calling the XF script from the command line is more difficult than it originally appeared to be. According to Remcom (the company that creates XF), you should be able to run an XF script from the command line with the option --execute-macro-script=/path/to/macro. However, this isn't supported (we think) by the version of XF that we have on OSC, and so they are looking in (talk to Carl abou this) to updating XF and how much that would cost. I'm not entirely sure that this is a solution either, because this requires calling the GUI for XF and I'm not sure that this is able to be done in a batch submission (clusters don't like using cluster computing with GUIs). Thus, it might be worthwhile to look into using NEC or HFSS, which I don't know anything about.
Have fun and let me know if you need any clarification or help,
Jordan |
/*
* File: Evolved_Dipole.cpp
* Author: Jordan Potter
*
* Created on July 14, 2017, 11:07 AM
*/
#include <cstdlib>
#include <iostream>
#include <cmath>
#include <fstream>
#include <string>
#include <sstream>
using namespace std;
const int NVAR=2; //number of variables (r and l)
const int NPOP=5; //size of population
double rand_normal(double mean, double stddev)
{//Box muller method to generate normal numbers
static double n2 = 0.0;
static int n2_cached = 0;
if (!n2_cached)
{
double x, y, r;
do
{
x = 2.0*rand()/RAND_MAX - 1;
y = 2.0*rand()/RAND_MAX - 1;
r = x*x + y*y;
}
while (r == 0.0 || r > 1.0);
{
double d = sqrt(-2.0*log(r)/r);
double n1 = x*d;
n2 = y*d;
double result = n1*stddev + mean;
n2_cached = 1;
return result;
}
}
else
{
n2_cached = 0;
return n2*stddev + mean;
}
}
void GenerateSample(double new_val[][NVAR], double mean_val[], double std_dev[],int n){
double u;
for(int i=0;i<NVAR;i++)
{
//Create a sample of mean + stdev*(gaussian sample) to use
u=rand_normal(0,.5);
new_val[n][i]=mean_val[i]+std_dev[i]*u;
//Include conditionals for bounds and constraints here
int errors=0;
if(new_val[n][i] <= 0)
{
i--; //If a parameter is less than 0, decrement so this value is created again
errors++;
}
if(errors >= 50) //If there are 50 errors something is going wrong, so stop
{
cout << "Error: Sample Generation Failed" << endl;
return;
}
}
}
void Simulation(double x[][NVAR], double y[NVAR], double mean[], double deviation[], double best[]){
//Write to CSV file. Might need to slightly adjust for generations after #1...
//For now, this will only write the last generation to file
ofstream handshake;
handshake.open("handshake.csv");
handshake << "C++ ESTRA -- Jordan Potter" << "\n";
handshake << "Mean Vector:" << "\n"; //Write the Current Mean Values to handshake.csv
for(int i=0;i<NVAR;i++)
{
if(i==(NVAR-1))
{
handshake << mean[i] << "\n";
}
else
{
handshake << mean[i] << ",";
}
}
handshake << "Deviation Vector:" << "\n"; //Write the Current Deviation Values to handshake.csv
for(int i=0;i<NVAR;i++)
{
if(i==(NVAR-1))
{
handshake << deviation[i] << "\n";
}
else
{
handshake << deviation[i] << ",";
}
}
handshake << "Best Antenna Vector:" << "\n"; //Write the Current Deviation Values to handshake.csv
for(int i=0;i<(NVAR+1);i++)
{
if(i==(NVAR))
{
handshake << best[i] << "\n";
}
else
{
handshake << best[i] << ",";
}
}
handshake << "Generation:" << "\n"; //Write the Generation Values to handshake.csv
for(int i=0;i<NPOP;i++)
{
for(int j=0;j<NVAR;j++)
{
if(j==(NVAR-1))
{
handshake << x[i][j] << "\n";
}
else
{
handshake << x[i][j] << ",";
}
}
}
handshake.close();
//Set Value in another file, so that controller.sh knows to start XF
}
void Read(double x[][NVAR], double y[NVAR], double mean[], double deviations[], double best[]){
//Import values from Handshake for mean/x population/current deviations
ifstream handshake;
handshake.open("handshake.csv");
string csvContent[NPOP+8]; //contain each line of csv
string strToDbl; //gets overwritten to contain each cell of a line. Then transferred to x,y,mean,deviations,best
for(int i=0;i<(NPOP+8);i++)
{
getline(handshake,csvContent[i]); //read in mean
if(i==2)
{
istringstream stream(csvContent[i]);
for(int j=0;j<NVAR;j++)
{
getline(stream,strToDbl,',');
mean[j] = atof(strToDbl.c_str());
}
}
else if(i==4) //read in deviations
{
istringstream stream(csvContent[i]);
for(int j=0;j<NVAR;j++)
{
getline(stream,strToDbl,',');
deviations[j] = atof(strToDbl.c_str());
}
}
else if(i==6) //read in the current best values
{
istringstream stream(csvContent[i]);
for(int j=0;j<NVAR+1;j++)
{
getline(stream,strToDbl,',');
best[j] = atof(strToDbl.c_str());
}
}
else if(i>=8) //read in the generation
{
istringstream stream(csvContent[i]);
for(int j=0;j<NVAR;j++)
{
getline(stream,strToDbl,',');
x[i-8][j] = atof(strToDbl.c_str());
}
}
}
//Import Values from Simulation from Handshook that XF (or eventually ARASim) will write
//May need to adjust if format of handshook changes
ifstream handshook;
handshook.open("handshook.csv");
string strToDbl2[NPOP+2];
for(int i=0;i<(NPOP+2);i++)
{
getline(handshook,strToDbl2[i]);
if(i>=2)
{
y[i-2] = atof(strToDbl2[i].c_str());
}
}
handshook.close();
}
int Mutate(double x[][NVAR],double out[], double best[], double mean[]){
//check to see if the output from x[i][] is larger than the output from x[m].
//i.e. output[i]>mean_output
int count=0;
for(int i=0;i<NPOP;i++)
{
if(out[i]>best[0])
{
best[0]=out[i];
for(int j=0;j<NVAR;j++)
{
best[j+1]=x[i][j];
mean[j]=x[i][j];
}
count++;
}
}
return count;
}
void CheckConvergence(int numSuccess, double p, double q, double d[]){
//If things aren't getting better than expand search radius.
//However, if things are getting better tighten search radius.
//Check this over k? values of x.
//P(output)>p then d=d/q
//if not then d=d*q
double P = (double) numSuccess/NPOP;
if(P >= p)
{
for(int i=0;i<NVAR;i++)
{
d[i]=d[i]/q;
}
}
else
{
for(int i=0;i<NVAR;i++)
{
d[i]=d[i]*q;
}
}
}
int main(int argc, char** argv) {
double m[NVAR]={.5,.5}; //mean values of r and l
double x[NPOP][NVAR]; //produced value of r and l
double d[NVAR] = {.25,.25}; //standard deviation vector
double best[NVAR+1]={-40,0,0}; //array to store param values and output score of top antenna {fit_score,r,l}
double output[NPOP]; //array to store scores of each generations
double q = .9; //factor to adjust search radius by
double p = .2; //ratio of samples that must be correct to adjust search radius
int numSuccess; //number of samples that improve on the best value
srand(time(0));
if(argc != 2)
cout << "Error: Usage. Specify --start or --cont" << endl;
else
{
if(string(argv[1]) == "--start")
{
for(int i=0;i<NPOP;i++)
{
GenerateSample(x,m,d,i); //Generate NPOP samples
}
Simulation(x, output, m, d, best); // Write current population to disk as well as mean,deviations and best population
}
else if(string(argv[1]) == "--cont")
{
//Read Function (read in population to x[],mean to m[],deviations to d[], best to best[] and scores to output[].
Read(x,output,m,d,best);
numSuccess = Mutate(x, output, best, m);
CheckConvergence(numSuccess, p, q, d);
//Check if d/d_o<lambda has been achieved. Consider adding this to check convergence function.
cout << "The best parameter values are:" << endl;
for(int i =0; i <NVAR; i++)
{
cout << m[i] << endl; //Consider writing to a log file
}
for(int i=0;i<NPOP;i++)
{
GenerateSample(x,m,d,i);
}
Simulation(x, output, m, d, best);
}
else
{
cout << "Error: Specify --start or --cont" << endl;
}
}
return 0;
}
|
var query = new ResultQuery();
///////////////////////Get Theta and Phi Gain///////////////
query.projectId = App.getActiveProject().getProjectDirectory();
query.simulationId = "000062";
query.runId = "Run0001";
query.sensorType = ResultQuery.FarZoneSensor;
query.sensorId = "Far Zone Sensor";
query.timeDependence = ResultQuery.SteadyState;
query.resultType = ResultQuery.Gain;
query.fieldScatter = ResultQuery.TotalField;
query.resultComponent = ResultQuery.Theta;
query.dataTransform = ResultQuery.NoTransform;
query.complexPart = ResultQuery.NotComplex;
query.surfaceInterpolationResolution = ResultQuery.NoInterpolation;
query.setDimensionRange( "Frequency", 0, -1 );
query.setDimensionRange( "Theta", 0, -1 );
query.setDimensionRange( "Phi", 0, -1 );
var thdata = new ResultDataSet( "" );
thdata.setQuery( query );
if( !thdata.isValid() ){
Output.println( "1getCurrentDataSet() : " +
thdata.getReasonWhyInvalid() );
}
query.resultComponent = ResultQuery.Phi;
var phdata = new ResultDataSet("");
phdata.setQuery(query);
if( !phdata.isValid() ){
Output.println( "2getCurrentDataSet() : " +
phdata.getReasonWhyInvalid() );
}
/////////////////Get Theta and Phi Phase///////////////////////////////////
query.resultType = ResultQuery.E;
query.fieldScatter = ResultQuery.TotalField;
query.resultComponent = ResultQuery.Theta;
query.dataTransform = ResultQuery.NoTransform;
query.complexPart = ResultQuery.Phase;
query.surfaceInterpolationResolution = ResultQuery.NoInterpolation;
query.setDimensionRange( "Frequency", 0, -1 );
query.setDimensionRange( "Theta", 0, -1 );
query.setDimensionRange( "Phi", 0, -1 );
var thphase = new ResultDataSet("");
thphase.setQuery(query);
if( !thphase.isValid() ){
Output.println( "3getCurrentDataSet() : " +
thphase.getReasonWhyInvalid() );
}
query.resultComponent = ResultQuery.Phi;
query.ComplexPart = ResultQuery.Phase;
var phphase = new ResultDataSet("");
phphase.setQuery(query);
if( !phphase.isValid() ){
Output.println( "4getCurrentDataSet() : " +
phphase.getReasonWhyInvalid() );
}
/////////////////Get Input Power///////////////////////////
query.sensorType = ResultQuery.System;
query.sensorId = "System";
query.timeDependence = ResultQuery.SteadyState;
query.resultType = ResultQuery.NetInputPower;
query.fieldScatter = ResultQuery.NoFieldScatter;
query.resultComponent = ResultQuery.Scalar;
query.dataTransform = ResultQuery.NoTransform;
query.complexPart = ResultQuery.NotComplex;
query.surfaceInterpolationResolution = ResultQuery.NoInterpolation;
query.clearDimensions();
query.setDimensionRange("Frequency",0,-1);
var inputpower = new ResultDataSet("");
inputpower.setQuery(query);
if( !inputpower.isValid() ){
Output.println( "5getCurrentDataSet() : " +
inputpower.getReasonWhyInvalid() );
}
FarZoneUtils.exportToUANFile(thdata,thphase,phdata,phphase,inputpower,"/users/PAS0654/osu8742/Evolved_Dipole/data/1.uan");
|