# Script to build an amorphous silica system

# Create an options dialog
void createDefaultDialog(Dialog ui)
{
        Widget w;

        # Set title and vertical automatic filling of main dialog layout
        ui.title = "Amorphous Silica Builder";

	# Basic controls
	ui.verticalFill = TRUE;
	ui.addIntegerSpin("nSiliconSpin", "Number of silicon atoms", 1, 1000000, 100, 1000);
	w = ui.addIntegerSpin("nOxygenSpin", "Number of oxygen atoms", 1, 2000000, 100, 1000);
	w.enabled = FALSE;
	w = ui.addCheck("keepRatio", "Fixed Si/O (1:2) Ratio", TRUE);
	w.onInteger(0, 1, "sendinteger", "nOxygenSpin", "disabled");
	ui.verticalFill = FALSE;
	ui.addDoubleSpin("densitySpin", "System Density", 0.001, 1000.0, 0.1, 2.0);
	ui.addCombo("unitsCombo", "", "g/cm**3,atoms/Ang**3", 1);
}

# Show (execute) the dialog to allow user interaction. The show() function returns 0 for 'canceled' or '1' for 'ok'
if (!showDefaultDialog()) error("Dialog canceled.\n");

# Get widget values
Dialog ui = defaultDialog();
int nSilicon, nOxygen;
nSilicon = ui.asInteger("nSiliconSpin");
if (ui.asInteger("keepRatio")) nOxygen = nSilicon * 2;
else nOxygen = ui.asInteger("nOxygenSpin");

# Determine system size, assuming cubic cell
double l, volume;
if (ui.asInteger("unitsCombo") == 1)
{
	# Density provided in g/cm**3
	printf("Density of %f g/cm**3 requested\n", ui.asDouble("densitySpin"));
	double mass = (nSilicon*32.06)+(nOxygen*15.9994);
	volume = (mass / ui.asDouble("densitySpin")) * 1.0E24 / AVOGADRO;
}
else
{
	# Density provided in atoms/Ang**3
	printf("Density of %f atoms/Ang**3 requested\n", ui.asDouble("densitySpin"));
	volume = (nSilicon+nOxygen) / ui.asDouble("densitySpin");
}

printf("Volume of cell = %f Ang**3\n", volume);
l = volume^(1.0/3.0);
printf("Length of cubic cell = %f Ang\n", l);

# Create temporary forcefield
Forcefield sio2ff = newFF("SiO2");
sio2ff.units = "kj";
sio2ff.addType(1,"Si","Si",Si, "");
sio2ff.addType(2,"O","O",O, "");
sio2ff.addInter("lj", 1, 0.0, 0.1, 2.0);
sio2ff.addInter("lj", 2, 0.0, 0.1, 2.0);
sio2ff.finalise();

# Create new models and associate ff
Model msi = newModel("Si");
newAtom(Si);
setupComponent("number", 1, nSilicon);
Model mo = newModel("O");
newAtom(O);
setupComponent("number", 1, nOxygen);
newModel("Amorphous Silica");
cell(l,l,l,90,90,90);

disorder("None", TRUE);

# Clean-up
deleteModel(msi);
deleteModel(mo);
deleteFF(sio2ff);
