rtang/SOLARIS_Web_Simulation
1
0
Fork 0
Code Issues Pull Requests Projects Releases Wiki Activity

added EnergyLevelsPlot

Browse Source
This commit is contained in:
Ryan@fsunuc 2023-07-27 00:07:05 -04:00
parent 014a13340c
commit 9871bdeeaa
6 changed files with 690 additions and 387 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

69
EnergyLevelsPlot.html Normal file
View File

@ -0,0 +1,69 @@
<!DOCTYPE html>
<html>
<head>
<title>Energy Levels Plot</title>
<link rel="icon" type="image/x-icon" href="logos/SOLARIS_favicon.png">
<script src="https://cdn.plot.ly/plotly-2.16.1.min.js"></script>
<meta name="viewport" content="width=device-width, height=device-height, initial-scale=1.0, user-scalable=no, user-scalable=0"/>
</head>
<style>
body{
font-family: Arial, Helvetica, sans-serif;
background : #6DB33E;
}
.column{
float : left;
width: 650px;
padding: 0px;
}
.row:after {
content: "";
display: table;
clear: both;
}
.plotStyle{
width:650px;
height:600px;
}
hr {
height:4px;
background-color:#F7CF3C;
border-style:none;
border-width:none;
}
</style>
<body>
<h1>Energy Levels Plot</h1>
<table>
<tr>
<td style="text-align: right;"> Isotopes Name </td>
<td><Input type="text" style="width:60px" value="18O" id="isotopes_name" enterkeyhint="done"/></td>
</tr>
<tr>
<td style="text-align: right;"> Max Ex </td>
<td><Input type="text" style="width:60px" value="6" id="maxEx" enterkeyhint="done"/></td>
<td> MeV</td>
</tr>
<!-- <tr>
<td style="text-align: right;"> PlotRange </td>
<td><Input type="text" style="width:60px" value="6" id="plotRange" enterkeyhint="done"/></td>
<td> fm</td>
</tr> -->
<tr>
<td></td>
<td><button onclick="PlotLevels()" style="width:65px">Plot</button></td>
</tr>
</table>
<div id="Plot_Levels" class="plotStyle"></div>
<p></p>
</body>
<script src="EnergyLevelsPlot.js"></script>
</html>

151
EnergyLevelsPlot.js Normal file
View File

@ -0,0 +1,151 @@
var energy = [];
var jpi = [];
var Name;
function GetData(){
Name = document.getElementById('isotopes_name').value;
let maxEx = parseFloat(document.getElementById('maxEx').value);
console.log(maxEx);
let str = 'get_nuclear_data.py?isotopes_name=' + Name + "&maxEx=" + maxEx;
let client = new XMLHttpRequest();
client.onreadystatechange = function() {
let haha = client.responseText.split('\n').slice(17);
jpi = [];
energy = [];
haha.forEach(line =>{
// console.log(line);
if( line.includes("<tr><td style=") && line.length != 0) {
jpi.push(line.substring(98).slice(0,-10).trim());
energy.push(parseFloat(line.substring(43,54).trim())/1000.);
// console.log(jpi[jpi.length - 1] + ", " + energy[energy.length-1]);
}
});
}
client.open('GET', str, false);
client.send();
}
function PlotLevels(){
GetData();
Plotly.purge("Plot_Levels");
if( energy.length == 0 ) return;
// console.log( Name + " | num. states : " + energy.length);
const plotWidth = 300;
const plotHeight = 600;
const yMin = -1;
const maxExExp = Math.max(...energy);
// console.log(maxExExp);
// let maxY = parseFloat(document.getElementById('plotRange').value);
const fig = {
data: [],
layout: {
plot_bgcolor: 'white',
width: plotWidth,
height: plotHeight,
margin: { l: 0, r: 0, t: 0, b: 0 },
showlegend: false,
xaxis: {
showline: false,
visible: false,
range: [-0.5, 3]
},
yaxis: {
range: [yMin, maxExExp + 1],
showline: false,
visible: false
},
annotations: []
}
};
const l = energy.length;
const fontSize = 14;
const fontSizeMeV = fontSize / plotHeight * (maxExExp + 1 - yMin);
let ypos = [];
for( let i = 0; i < energy.length; i++) ypos.push(energy[i]);
let noOverlap = false;
let loop = 0;
while (!noOverlap && loop < 2 * l) {
for (let i = 1; i <= l; i++) {
const diff = ypos[i] - ypos[i - 1];
if (diff < fontSizeMeV) {
ypos[i - 1] += (diff - fontSizeMeV) / 2;
ypos[i] += (fontSizeMeV - diff) / 2;
if (ypos[i - 1] < yMin + fontSizeMeV / 2) {
ypos[i - 1] = yMin + fontSizeMeV / 2;
ypos[i] = ypos[i - 1] + fontSizeMeV;
}
}
}
let count = 0;
for (let i = 1; i <= l; i++) {
const diff = ypos[i] - ypos[i - 1];
if (diff > fontSizeMeV) {
count++;
}
}
if (count === l) {
noOverlap = true;
}
loop++;
}
for (let i = 0; i < l; i++) {
fig.data.push({
x: [0, 1],
y: [energy[i], energy[i]],
mode: 'lines',
line: { color: 'black', width: 1 }
});
fig.data.push({
x: [1.03, 1.1, 1.19],
y: [energy[i], ypos[i], ypos[i]],
mode: 'lines',
line: { color: 'gray', width: 1 }
});
// console.log(energy[i]+ ", " + ypos[i]);
fig.layout.annotations.push({
x: 1.2,
y: ypos[i],
text: `${energy[i].toFixed(3)}, ${jpi[i]}`,
xanchor: 'left',
font: { size: fontSize },
showarrow: false
});
}
// let NameYPos = (parseFloat(maxEx) + 2*fontSizeMeV);
// console.log(NameYPos);
fig.layout.annotations.push({
x: 0.5,
y: (parseFloat(maxEx) + 0.5),
text: Name,
font: { size: 2 * fontSize },
showarrow: false
});
// Create the plot
Plotly.newPlot('Plot_Levels', fig.data, fig.layout);
}

6
get_nuclear_data.py
View File

@ -16,8 +16,8 @@ def lc_read_csv(url):
req.add_header('User-Agent', 'Mozilla/5.0 (X11; Ubuntu; Linux x86_64; rv:77.0) Gecko/20100101 Firefox/77.0')
return pd.read_csv(urllib.request.urlopen(req))
mp = 938.27208816; #MeV/c^2
mn = 939.56542052;
mp = 938.27208816 #MeV/c^2
mn = 939.56542052
#===================================================
import cgi, cgitb
@ -31,7 +31,7 @@ if maxEx == "can be omitted" :
maxEx = -1
query = livechart + "fields=ground_states&nuclides=" + AZ
temp = lc_read_csv(query);
temp = lc_read_csv(query)
print( "Content-type:text/html\r\n\r\n")
print("<html>")

22
index.html
View File

@ -7,6 +7,7 @@
</head>
<meta charset="utf-8">
<meta name="viewport" content="width=device-width, initial-scale=1">
<script src="https://cdn.plot.ly/plotly-2.16.1.min.js"></script>
<style>
:root{
@ -124,7 +125,7 @@
<td style="text-align:right"><a href="instruction.html" target="uploaded">Intructions & Credits</a></td>
</tr>
<tr>
<td style="text-align:right"><a href="test.html" target="_blank">Woods-Saxon (const.)</a></td>
<td style="text-align:right"><a href="test.html" target="_blank">Woods-Saxon</a></td>
</tr>
<tr>
<td style="text-align:right"><a href="nuclearChart.html" target="_blank">Nuclides Chart (const.)</a></td>
@ -137,21 +138,25 @@
<table class="center">
<tr>
<td style="text-align:right">Isotopes Name:</td>
<td><input type = "text" name = "isotopes_name" size="13" value="24F" enterkeyhint="done"/></td>
<td><input type = "text" name = "isotopes_name" id="isotopes_name" size="13" value="24F" enterkeyhint="done"/></td>
</tr>
<tr>
<td style="text-align:right">Max Ex [MeV]:</td>
<td><input type = "text" name = "maxEx" size="13" value="can be omitted" enterkeyhint="done"/></td>
<td><input type = "text" name = "maxEx" id="maxEx" size="13" value="can be omitted" enterkeyhint="done"/></td>
</tr>
<tr>
<!-- <tr>
<td></td>
<td>
<input type = "submit" value = "Get Isotopes Data" />
<input type = "submit" value = "Get Isotopes Data"/>
</td>
</tr>
</tr> -->
<tr>
<td></td>
<td><button onclick="PlotLevels()" style="width: 130px;">Get Isotopes Data</button></td>
</tr>
</table>
</form>
<div id="Plot_Levels" class="plotStyle"></div>
<iframe name="NuclearData" style="border:none;width:100%;height:100%" ></iframe>
</nav>
@ -162,6 +167,9 @@
</article>
</section>
<script src="EnergyLevelsPlot.js"></script>
<script>
/*

4
simpleSim.html
View File

@ -281,11 +281,11 @@ Max Ex: <input type="text" id="maxEx" size="5" value="5"/>MeV
<td><input type = "file" name = "filename4a" /> </td>
</tr>
<tr>
<td style="text-align:right" width="200">Plot Config File * </th>
<td style="text-align:right" width="200">Plot Config File # </th>
<td><input type = "file" name = "filename5" /> </td>
</tr>
<tr>
<td>^ can be alone <br>* can be omitted</th>
<td>^ can be alone <br># can be omitted</th>
<td><input type = "submit" value = "Upload & Run Simulation" style="height:50px; width:200px" formtarget="_blank"/> </td>
</tr>
</table>

825
test.cpp
View File

@ -18,6 +18,7 @@
#include "implot.h"
#include <iostream>
#include <cstring>
#include <vector>
#include "WS.h"
@ -92,414 +93,488 @@ void loop(){
ImGui::End();
}
// Woods-Saxon window
{
ImGui::SetNextWindowSize(ImVec2(600, 1000), ImGuiCond_FirstUseEver);
ImGui::SetNextWindowPos(ImVec2(50, 50), ImGuiCond_FirstUseEver);
static float V0 = -45, r0 = 1.25, R0 = 3.5, a0 = 0.6;
static float VSO = 28, rSO = 1.25, RSO = 3.5, aSO = 0.6;
static int Z = 0, nStep = 300;
static float Rc = 3.5, rc = 1.25, dr = 0.1;
static int A = 20;
static float mass = 939.565;
ImGui::Begin("Woods-Saxon Calculation");
static int e = 0;
if( ImGui::RadioButton("Arbitary", &e, 0) ){
A = 1;
}
ImGui::SameLine();
if( ImGui::RadioButton("Isotope", &e, 1) ){
rc = 1.25;
r0 = 1.25;
rSO = 1.25;
}
if( e == 0){
rc = Rc / pow(A, 1./3);
r0 = R0 / pow(A, 1./3);
rSO = RSO / pow(A, 1./3);
}
if( e == 1){
Rc = rc * pow(A, 1./3);
R0 = r0 * pow(A, 1./3);
RSO = rSO * pow(A, 1./3);
}
ImGui::Separator();
ImGui::PushItemWidth(100);
static int nucleon = 0;
ImGui::RadioButton("Neutron", &nucleon, 0); ImGui::SameLine();
ImGui::RadioButton("Proton", &nucleon, 1);
ImGui::BeginDisabled(e==0); ImGui::DragInt("A ", &A, 1, 1, 200); ImGui::EndDisabled();
ImGui::SameLine();
ImGui::DragInt("Z", &Z, 1, 0, 100);
ImGui::BeginDisabled(e==1); ImGui::DragFloat("Rc [fm] ", &Rc, 0.01, 1, 10); ImGui::EndDisabled();
ImGui::SameLine();
ImGui::BeginDisabled(e==0);
ImGui::DragFloat("rc [fm] ", &rc, 0.01, 1, 10);
ImGui::EndDisabled();
ImGui::Text("Eff. mass : %.3f ", mass);
ImGui::Separator();
ImGui::DragFloat("V0 [MeV]", &V0, 0.1, -100, 0);
ImGui::BeginDisabled(e==1); ImGui::DragFloat("R0 [fm] ", &R0, 0.01, 1, 10); ImGui::EndDisabled();
ImGui::SameLine();
ImGui::BeginDisabled(e==0); ImGui::DragFloat("r0 [fm] ", &r0, 0.01, 1, 10); ImGui::EndDisabled();
ImGui::DragFloat("a0 [fm]", &a0, 0.01, 0.1, 2);
ImGui::DragFloat("VSO [MeV]", &VSO, 0.1, 0, 40);
ImGui::BeginDisabled(e==1); ImGui::DragFloat("RSO [fm] ", &RSO, 0.01, 1, 10); ImGui::EndDisabled();
ImGui::SameLine();
ImGui::BeginDisabled(e==0); ImGui::DragFloat("rSO [fm] ", &rSO, 0.01, 1, 10); ImGui::EndDisabled();
ImGui::DragFloat("aSO [fm]", &aSO, 0.01, 0.1, 2);
ImGui::Separator();
ImGui::DragInt("nStep", &nStep, 50, 100, 400);
ImGui::SameLine();
ImGui::DragFloat("dr [fm]", &dr, 0.001, 0.01, 0.1);
ImGui::PushStyleVar(ImGuiStyleVar_FramePadding, ImVec2(30, 10));
if( ImGui::Button("Cal WS Levels") ){
WoodsSaxon ws;
if( nucleon == 0) {
ws.IsNeutron();
}else{
ws.IsProton();
}
if( e == 0) {
ws.SetNucleus(1,Z);
}else{
ws.SetNucleus(A,Z);
mass = 931.5 * (ws.GetA() )/(1.008664 + ws.GetA());
ws.SetMass(mass);
}
if( Z > 0 ) ws.SetRc(Rc);
ws.SetV0( V0 );
ws.SetR0( R0 );
ws.Seta0( a0 );
ImGui::SetNextWindowSize(ImVec2(600, 1000), ImGuiCond_FirstUseEver);
ImGui::SetNextWindowPos(ImVec2(50, 50), ImGuiCond_FirstUseEver);
ws.SetVSO( VSO );
ws.SetRSO( RSO );
ws.SetaSO( aSO );
ws.SetRange2(0.0001, dr, nStep);
ws.CalWSEnergies(false, 7, 100, 1e-7, 50, 0.2, false);
float dx = nStep * dr / nPt;
wfr.clear();
for( int i = 0; i < nPt; i ++){
float r = i * dx;
xValues[i] = r;
WSCValues[i] = V0/(1 + exp(( r - R0)/a0));
WSSOValues[i] = VSO * exp((r-RSO)/aSO) / pow(1+exp((r-RSO)/aSO), 2) / aSO/ r ;
wfr.push_back(dr*i);
}
wf.clear();
for( int i = 0; i < ws.orbString.size(); i++){
wf.push_back(ws.CalWaveFunction(i, abs(V0)/2, ws.energy[i]));
}
selected = -1;
energies.clear();
orbString.clear();
energies = ws.energy;
orbString = ws.orbString;
//ws.PrintEnergyLevels();
}
ImGui::PopStyleVar();
ImGui::SetNextItemOpen(true, ImGuiCond_Always);
if( ImGui::TreeNode("WS Energy Levels:") ){
for( int i = 0; i < energies.size(); i++){
char buf[32];
sprintf(buf, "%24.12f MeV %s", energies[i],orbString[i].c_str());
if( ImGui::Selectable(buf, selected == i) ) selected = i;
}
ImGui::TreePop();
}
if( ImPlot::BeginPlot("Plot", ImVec2(ImGui::GetWindowWidth()-20, 400)) ){
ImPlot::SetupLegend(ImPlotLocation_SouthEast);
ImPlot::SetupAxes("r [fm]"," Energy [MeV]");
ImPlot::SetupAxesLimits(0, 10, floor(V0*1.1), 10);
ImPlot::SetupAxisLimitsConstraints(ImAxis_X1, 0, 30);
ImPlot::SetupAxisLimitsConstraints(ImAxis_Y1, floor(V0*1.1), abs(floor(V0*1.1)));
ImPlot::PlotLine("Central", xValues, WSCValues, nPt);
ImPlot::PlotLine("S-O", xValues, WSSOValues, nPt);
if( selected >= 0 ){
for( int i = 0; i < nPt; i ++) Energy[i] = energies[selected];
ImPlot::PlotLine(orbString[selected].c_str(), xValues, Energy, nPt);
const double * haha = wf[selected].data();
const double * kaka = wfr.data();
ImPlot::PlotLine("WF", kaka, haha, wfr.size());
}
ImPlot::EndPlot();
}
ImGui::End();
}
//*============= WS Trend
{
ImGui::SetNextWindowSize(ImVec2(600, 800), ImGuiCond_FirstUseEver);
ImGui::SetNextWindowPos(ImVec2(680, 50), ImGuiCond_FirstUseEver);
static float V0 = -45, r0 = 1.25, a0 = 0.6;
static float VSO = 28, rSO = 1.25, aSO = 0.6;
static int Z = 0, N = 1, nStep = 300;
static float rc = 1.25, dr = 0.1;
static int A[2] = {10, 20};
static int NRange[2] = {1, 20};
static int ZRange[2] = {1, 20};
static float mass = 939.565;
static float kappa = -1;
ImGui::Begin("Woods-Saxon Calculation");
if( ImGui::BeginTabBar("Woods-Saxon", ImGuiTabBarFlags_None)){
if (ImGui::BeginTabItem("Cal.")){
ImGui::Begin("Woods-Saxon Calculation (Range)");
static float V0 = -45, r0 = 1.25, R0 = 3.5, a0 = 0.6;
static float VSO = 28, rSO = 1.25, RSO = 3.5, aSO = 0.6;
static int Z = 0, nStep = 300;
static float Rc = 3.5, rc = 1.25, dr = 0.1;
static int A = 20;
static float mass = 939.565;
static int e = 0;
ImGui::RadioButton("Stablility", &e, 0); ImGui::SameLine();
ImGui::RadioButton("var. A", &e, 1); ImGui::SameLine();
ImGui::RadioButton("var. N", &e, 2); ImGui::SameLine();
ImGui::RadioButton("var. Z", &e, 3);
ImGui::Separator();
static int e = 0;
if( ImGui::RadioButton("Arbitary", &e, 0) ){
A = 1;
}
ImGui::SameLine();
if( ImGui::RadioButton("Isotope", &e, 1) ){
rc = 1.25;
r0 = 1.25;
rSO = 1.25;
}
ImGui::PushItemWidth(150);
static int nucleon = 0; //neutron or proton
if( e != 1 ){
ImGui::RadioButton("Neutron", &nucleon, 0); ImGui::SameLine();
ImGui::RadioButton("Proton", &nucleon, 1);
}else{
ImGui::Text("Neutron only");
}
if( e == 0 || e == 1 ){
ImGui::DragInt2("A range", A, 1, 1, 300);
}
if( e == 2){
ImGui::DragInt("Z", &Z, 0, 0, 100);
ImGui::DragInt2("N range", NRange, 1, 1, 300);
}
if( e == 3){
ImGui::DragInt("N", &N, 1, 1, 200);
ImGui::DragInt2("Z range", ZRange, 1, 0, 300);
}
if( e != 1){
ImGui::DragFloat("kappa", &kappa, 0.1, -10, 10);
}
ImGui::Separator();
ImGui::DragFloat("V0 [MeV]", &V0, 0.1, -100, 0);
ImGui::DragFloat("r0 [fm] ", &r0, 0.01, 1, 10);
ImGui::DragFloat("a0 [fm]", &a0, 0.01, 0.1, 2);
ImGui::DragFloat("VSO [MeV]", &VSO, 0.1, 0, 40);
ImGui::DragFloat("rSO [fm] ", &rSO, 0.01, 1, 10);
ImGui::DragFloat("aSO [fm]", &aSO, 0.01, 0.1, 2);
ImGui::Separator();
ImGui::DragInt("nStep", &nStep, 50, 100, 400);
ImGui::SameLine();
ImGui::DragFloat("dr [fm]", &dr, 0.001, 0.01, 0.1);
int ARange[2] ;
int fixedNZ;
if( e == 0 || e == 1 ) {ARange[0] = A[0]; ARange[1] = A[1]; fixedNZ = 0;}
if( e == 2 ) {ARange[0] = NRange[0] + Z; ARange[1] = NRange[1] + Z; fixedNZ = Z;}
if( e == 3 ) {ARange[0] = ZRange[0] + N; ARange[1] = ZRange[1] + N; fixedNZ = N;}
int maxA = ARange[1];
if( maxA < A[0]) {
maxA = ARange[0];
ARange[0] = ARange[1];
ARange[1] = maxA;
}
orbitalStr.clear();
orbitalStr.push_back("0s1/2"); ///2
orbitalStr.push_back("0p3/2");
orbitalStr.push_back("0p1/2"); ///8
if( maxA > 8 ) {
orbitalStr.push_back("0d5/2");
orbitalStr.push_back("1s1/2");
orbitalStr.push_back("0d3/2"); ///20
}
if( maxA > 20 ){
orbitalStr.push_back("0f7/2"); ///28
orbitalStr.push_back("1p3/2");
orbitalStr.push_back("0f5/2");
}
if( maxA > 40 ){
orbitalStr.push_back("1p1/2"); ///40
orbitalStr.push_back("0g9/2"); ///50
orbitalStr.push_back("0g7/2");
orbitalStr.push_back("1d5/2");
orbitalStr.push_back("0h11/2");
orbitalStr.push_back("1d3/2");
}
if( maxA > 82 ){
orbitalStr.push_back("2s1/2"); ///82
orbitalStr.push_back("0h9/2");
orbitalStr.push_back("1f7/2");
orbitalStr.push_back("0i13/2");
orbitalStr.push_back("2p3/2");
orbitalStr.push_back("1f5/2");
}
if( maxA > 126 ){
orbitalStr.push_back("2p1/2"); ///126
orbitalStr.push_back("1g9/2");
orbitalStr.push_back("0i11/2");
orbitalStr.push_back("0j15/2");
orbitalStr.push_back("0j15/2");
orbitalStr.push_back("2d5/2");
orbitalStr.push_back("3s1/2");
orbitalStr.push_back("2d3/2");
orbitalStr.push_back("1g7/2"); ///184
}
static float progress = 0;
static bool cal = false;
static bool finsihed = false;
static int ANZ = ARange[0];
if( ImGui::Button("Cal WS Trend") ){
cal = true;
finsihed = false;
AList.clear();
enLevelList.clear();
ANZ = ARange[0];
}
if(cal) {
WoodsSaxon ws;
ws.Setr0(r0);
ws.SetrSO(rSO);
ws.Setrc(rc);
ws.SetVSO(VSO);
ws.SetRange2(0.001, dr, nStep);
if( nucleon == 0 ) {
ws.IsNeutron();
}else{
ws.IsProton();
}
if( ANZ <= ARange[1] ){
AList.push_back(ANZ);
if( e == 0){
ws.SetNucleus ( ANZ, 0.008 + 0.495 * ANZ - 0.00076 * ANZ *ANZ + 1.4e-6 * ANZ * ANZ * ANZ); // aproximate stable vallage
ws.CalRadius();
ws.SetV0 ( V0 * ( 1 + kappa * ( 2.0* ws.GetN() - ws.GetA() ) / ws.GetA() ) );
rc = Rc / pow(A, 1./3);
r0 = R0 / pow(A, 1./3);
rSO = RSO / pow(A, 1./3);
}
if( e == 1){
ws.SetNucleus ( ANZ, 0);
ws.CalRadius();
ws.SetV0( V0 );
Rc = rc * pow(A, 1./3);
R0 = r0 * pow(A, 1./3);
RSO = rSO * pow(A, 1./3);
}
ImGui::Separator();
if( e == 2 ){ /// change in N, ANZ = N
ws.SetNucleus( ANZ , fixedNZ );
ws.CalRadius();
ImGui::PushItemWidth(100);
static int nucleon = 0;
ImGui::RadioButton("Neutron", &nucleon, 0); ImGui::SameLine();
ImGui::RadioButton("Proton", &nucleon, 1);
ImGui::BeginDisabled(e==0); ImGui::DragInt("A ", &A, 1, 1, 200); ImGui::EndDisabled();
ImGui::SameLine();
ImGui::DragInt("Z", &Z, 1, 0, 100);
ws.SetV0 ( V0 * ( 1 + kappa * ( 2.0* ws.GetN() - ws.GetA() ) / ws.GetA() ) );
}
if( e == 3 ){ /// change in Z, ANZ = Z
ws.SetNucleus( ANZ , ANZ - fixedNZ );
ws.CalRadius();
ws.SetV0 ( V0 * ( 1 + kappa * ( 2.0* ws.GetN() - ws.GetA() ) / ws.GetA() ) );
}
ImGui::BeginDisabled(e==1); ImGui::DragFloat("Rc [fm] ", &Rc, 0.01, 1, 10); ImGui::EndDisabled();
ImGui::SameLine();
ImGui::BeginDisabled(e==0);
ImGui::DragFloat("rc [fm] ", &rc, 0.01, 1, 10);
ImGui::EndDisabled();
double reducedMass = 931.5 * (ws.GetA() )/(1.008664 + ws.GetA());
ws.SetMass(reducedMass);
ImGui::Text("Eff. mass : %.3f ", mass);
ws.ClearVector();
ws.CalWSEnergies(false, 7, 500, 1e-7, 400, 0.2, false);
ImGui::Separator();
//ws.PrintWSParas();
// ws.PrintEnergyLevels();
ImGui::DragFloat("V0 [MeV]", &V0, 0.1, -100, 0);
vector<double> exList;
for( int i = 0; i < orbitalStr.size(); i++){
bool isMatched = false;
for( int j = 0; j < (ws.energy).size(); j++){
if( orbitalStr[i] == ws.orbString[j] ) {
exList.push_back(ws.energy[j]);
isMatched = true;
}
ImGui::BeginDisabled(e==1); ImGui::DragFloat("R0 [fm] ", &R0, 0.01, 1, 10); ImGui::EndDisabled();
ImGui::SameLine();
ImGui::BeginDisabled(e==0); ImGui::DragFloat("r0 [fm] ", &r0, 0.01, 1, 10); ImGui::EndDisabled();
ImGui::DragFloat("a0 [fm]", &a0, 0.01, 0.1, 2);
ImGui::DragFloat("VSO [MeV]", &VSO, 0.1, 0, 40);
ImGui::BeginDisabled(e==1); ImGui::DragFloat("RSO [fm] ", &RSO, 0.01, 1, 10); ImGui::EndDisabled();
ImGui::SameLine();
ImGui::BeginDisabled(e==0); ImGui::DragFloat("rSO [fm] ", &rSO, 0.01, 1, 10); ImGui::EndDisabled();
ImGui::DragFloat("aSO [fm]", &aSO, 0.01, 0.1, 2);
ImGui::Separator();
ImGui::DragInt("nStep", &nStep, 50, 100, 400);
ImGui::SameLine();
ImGui::DragFloat("dr [fm]", &dr, 0.001, 0.01, 0.1);
ImGui::PushStyleVar(ImGuiStyleVar_FramePadding, ImVec2(30, 10));
if( ImGui::Button("Cal WS Levels") ){
WoodsSaxon ws;
if( nucleon == 0) {
ws.IsNeutron();
}else{
ws.IsProton();
}
if( isMatched == false ) exList.push_back(0);
if( e == 0) {
ws.SetNucleus(1,Z);
}else{
ws.SetNucleus(A,Z);
mass = 931.5 * (ws.GetA() )/(1.008664 + ws.GetA());
ws.SetMass(mass);
}
if( Z > 0 ) ws.SetRc(Rc);
ws.SetV0( V0 );
ws.SetR0( R0 );
ws.Seta0( a0 );
ws.SetVSO( VSO );
ws.SetRSO( RSO );
ws.SetaSO( aSO );
ws.SetRange2(0.0001, dr, nStep);
ws.CalWSEnergies(false, 7, 100, 1e-7, 50, 0.2, false);
float dx = nStep * dr / nPt;
wfr.clear();
for( int i = 0; i < nPt; i ++){
float r = i * dx;
xValues[i] = r;
WSCValues[i] = V0/(1 + exp(( r - R0)/a0));
WSSOValues[i] = VSO * exp((r-RSO)/aSO) / pow(1+exp((r-RSO)/aSO), 2) / aSO/ r ;
wfr.push_back(dr*i);
}
wf.clear();
for( int i = 0; i < ws.orbString.size(); i++){
wf.push_back(ws.CalWaveFunction(i, abs(V0)/2, ws.energy[i]));
}
selected = -1;
energies.clear();
orbString.clear();
energies = ws.energy;
orbString = ws.orbString;
//ws.PrintEnergyLevels();
}
enLevelList.push_back(exList);
ImGui::PopStyleVar();
progress = (ANZ - ARange[0]) * 1.0 / (ARange[1] - ARange[0]);
ImGui::SetNextItemOpen(true, ImGuiCond_Always);
if( ImGui::TreeNode("WS Energy Levels:") ){
for( int i = 0; i < energies.size(); i++){
char buf[32];
sprintf(buf, "%24.12f MeV %s", energies[i],orbString[i].c_str());
if( ImGui::Selectable(buf, selected == i) ) selected = i;
}
ImGui::TreePop();
}
ANZ ++;
}else{
cal = false;
finsihed = true;
if( ImPlot::BeginPlot("Plot", ImVec2(ImGui::GetWindowWidth()-20, 400)) ){
ImPlot::SetupLegend(ImPlotLocation_SouthEast);
ImPlot::SetupAxes("r [fm]"," Energy [MeV]");
ImPlot::SetupAxesLimits(0, 10, floor(V0*1.1), 10);
ImPlot::SetupAxisLimitsConstraints(ImAxis_X1, 0, 30);
ImPlot::SetupAxisLimitsConstraints(ImAxis_Y1, floor(V0*1.1), abs(floor(V0*1.1)));
ImPlot::PlotLine("Central", xValues, WSCValues, nPt);
ImPlot::PlotLine("S-O", xValues, WSSOValues, nPt);
if( selected >= 0 ){
for( int i = 0; i < nPt; i ++) Energy[i] = energies[selected];
ImPlot::PlotLine(orbString[selected].c_str(), xValues, Energy, nPt);
const double * haha = wf[selected].data();
const double * kaka = wfr.data();
ImPlot::PlotLine("WF", kaka, haha, wfr.size());
}
ImPlot::EndPlot();
}
ImGui::EndTabItem();
}
if( ImGui::BeginTabItem("Range")){
static float V0 = -45, r0 = 1.25, a0 = 0.6;
static float VSO = 28, rSO = 1.25, aSO = 0.6;
static int Z = 0, N = 1, nStep = 300;
static float rc = 1.25, dr = 0.1;
static int A[2] = {10, 20};
static int NRange[2] = {1, 20};
static int ZRange[2] = {1, 20};
static float mass = 939.565;
static float kappa = -1;
static int e = 0;
ImGui::RadioButton("Stablility", &e, 0); ImGui::SameLine();
ImGui::RadioButton("var. A", &e, 1); ImGui::SameLine();
ImGui::RadioButton("var. N", &e, 2); ImGui::SameLine();
ImGui::RadioButton("var. Z", &e, 3);
ImGui::Separator();
ImGui::PushItemWidth(150);
static int nucleon = 0; //neutron or proton
if( e != 1 ){
ImGui::RadioButton("Neutron", &nucleon, 0); ImGui::SameLine();
ImGui::RadioButton("Proton", &nucleon, 1);
}else{
ImGui::Text("Neutron only");
}
if( e == 0 || e == 1 ){
ImGui::DragInt2("A range", A, 1, 1, 300);
}
if( e == 2){
ImGui::DragInt("Z", &Z, 0, 0, 100);
ImGui::DragInt2("N range", NRange, 1, 1, 300);
}
if( e == 3){
ImGui::DragInt("N", &N, 1, 1, 200);
ImGui::DragInt2("Z range", ZRange, 1, 0, 300);
}
if( e != 1){
ImGui::DragFloat("kappa", &kappa, 0.1, -10, 10);
}
ImGui::Separator();
ImGui::DragFloat("V0 [MeV]", &V0, 0.1, -100, 0);
ImGui::DragFloat("r0 [fm] ", &r0, 0.01, 1, 10);
ImGui::DragFloat("a0 [fm]", &a0, 0.01, 0.1, 2);
ImGui::DragFloat("VSO [MeV]", &VSO, 0.1, 0, 40);
ImGui::DragFloat("rSO [fm] ", &rSO, 0.01, 1, 10);
ImGui::DragFloat("aSO [fm]", &aSO, 0.01, 0.1, 2);
ImGui::Separator();
ImGui::DragInt("nStep", &nStep, 50, 100, 400);
ImGui::SameLine();
ImGui::DragFloat("dr [fm]", &dr, 0.001, 0.01, 0.1);
int ARange[2] ;
int fixedNZ;
if( e == 0 || e == 1 ) {ARange[0] = A[0]; ARange[1] = A[1]; fixedNZ = 0;}
if( e == 2 ) {ARange[0] = NRange[0] + Z; ARange[1] = NRange[1] + Z; fixedNZ = Z;}
if( e == 3 ) {ARange[0] = ZRange[0] + N; ARange[1] = ZRange[1] + N; fixedNZ = N;}
int maxA = ARange[1];
if( maxA < A[0]) {
maxA = ARange[0];
ARange[0] = ARange[1];
ARange[1] = maxA;
}
orbitalStr.clear();
{
orbitalStr.push_back("0s1/2"); ///2
orbitalStr.push_back("0p3/2");
orbitalStr.push_back("0p1/2"); ///8
if( maxA > 8 ) {
orbitalStr.push_back("0d5/2");
orbitalStr.push_back("1s1/2");
orbitalStr.push_back("0d3/2"); ///20
}
if( maxA > 20 ){
orbitalStr.push_back("0f7/2"); ///28
orbitalStr.push_back("1p3/2");
orbitalStr.push_back("0f5/2");
}
if( maxA > 40 ){
orbitalStr.push_back("1p1/2"); ///40
orbitalStr.push_back("0g9/2"); ///50
orbitalStr.push_back("0g7/2");
orbitalStr.push_back("1d5/2");
orbitalStr.push_back("0h11/2");
orbitalStr.push_back("1d3/2");
}
if( maxA > 82 ){
orbitalStr.push_back("2s1/2"); ///82
orbitalStr.push_back("0h9/2");
orbitalStr.push_back("1f7/2");
orbitalStr.push_back("0i13/2");
orbitalStr.push_back("2p3/2");
orbitalStr.push_back("1f5/2");
}
if( maxA > 126 ){
orbitalStr.push_back("2p1/2"); ///126
orbitalStr.push_back("1g9/2");
orbitalStr.push_back("0i11/2");
orbitalStr.push_back("0j15/2");
orbitalStr.push_back("0j15/2");
orbitalStr.push_back("2d5/2");
orbitalStr.push_back("3s1/2");
orbitalStr.push_back("2d3/2");
orbitalStr.push_back("1g7/2"); ///184
}
}
static float progress = 0;
static bool cal = false;
static bool finsihed = false;
static int ANZ = ARange[0];
if( ImGui::Button("Cal WS Trend") ){
cal = true;
finsihed = false;
AList.clear();
enLevelList.clear();
ANZ = ARange[0];
}
if(cal) {
WoodsSaxon ws;
ws.Setr0(r0);
ws.SetrSO(rSO);
ws.Setrc(rc);
ws.SetVSO(VSO);
ws.SetRange2(0.001, dr, nStep);
if( nucleon == 0 ) {
ws.IsNeutron();
}else{
ws.IsProton();
}
if( ANZ <= ARange[1] ){
AList.push_back(ANZ);
if( e == 0){
ws.SetNucleus ( ANZ, 0.008 + 0.495 * ANZ - 0.00076 * ANZ *ANZ + 1.4e-6 * ANZ * ANZ * ANZ); // aproximate stable vallage
ws.CalRadius();
ws.SetV0 ( V0 * ( 1 + kappa * ( 2.0* ws.GetN() - ws.GetA() ) / ws.GetA() ) );
}
if( e == 1){
ws.SetNucleus ( ANZ, 0);
ws.CalRadius();
ws.SetV0( V0 );
}
if( e == 2 ){ /// change in N, ANZ = N
ws.SetNucleus( ANZ , fixedNZ );
ws.CalRadius();
ws.SetV0 ( V0 * ( 1 + kappa * ( 2.0* ws.GetN() - ws.GetA() ) / ws.GetA() ) );
}
if( e == 3 ){ /// change in Z, ANZ = Z
ws.SetNucleus( ANZ , ANZ - fixedNZ );
ws.CalRadius();
ws.SetV0 ( V0 * ( 1 + kappa * ( 2.0* ws.GetN() - ws.GetA() ) / ws.GetA() ) );
}
double reducedMass = 931.5 * (ws.GetA() )/(1.008664 + ws.GetA());
ws.SetMass(reducedMass);
ws.ClearVector();
ws.CalWSEnergies(false, 7, 500, 1e-7, 400, 0.2, false);
//ws.PrintWSParas();
// ws.PrintEnergyLevels();
vector<double> exList;
for( int i = 0; i < orbitalStr.size(); i++){
bool isMatched = false;
for( int j = 0; j < (ws.energy).size(); j++){
if( orbitalStr[i] == ws.orbString[j] ) {
exList.push_back(ws.energy[j]);
isMatched = true;
}
}
if( isMatched == false ) exList.push_back(0);
}
enLevelList.push_back(exList);
progress = (ANZ - ARange[0]) * 1.0 / (ARange[1] - ARange[0]);
ANZ ++;
}else{
cal = false;
finsihed = true;
}
}
ImGui::SameLine(); ImGui::ProgressBar(progress);
if( finsihed ){
if( ImPlot::BeginPlot("Plot", ImVec2(ImGui::GetWindowWidth()-20, 400)) ){
ImPlot::SetupLegend(ImPlotLocation_East, ImPlotLegendFlags_Outside);
ImPlot::SetupAxes("A"," Energy [MeV]");
ImPlot::SetupAxesLimits(ARange[0], ARange[1], -50, 0);
// ImPlot::SetupAxisLimitsConstraints(ImAxis_X1, ARange[0], ARange[1]);
// ImPlot::SetupAxisLimitsConstraints(ImAxis_Y1, floor(V0)*2, 0);
const double * kaka = AList.data();
for( int i = 0; i < orbitalStr.size(); i++){
vector<double> jaja;
for( int j = 0; j < AList.size(); j++) jaja.push_back(enLevelList[j][i]);
const double * haha = jaja.data();
ImPlot::PlotLine(orbitalStr[i].c_str(), kaka, haha, AList.size());
}
ImPlot::EndPlot();
}
static bool saveDataClicked = false;
if( ImGui::Button("Download data") ){
// FILE * file_out;
// file_out = fopen("", "w+");
saveDataClicked = true;
}
if( saveDataClicked ) {
ImGui::Text("to be impletmented....");
}
}
ImGui::EndTabItem();
}
if( ImGui::BeginTabItem("Fit")){
ImGui::Text("to be impletment....");
ImGui::EndTabItem();
}
ImGui::EndTabBar();
}
ImGui::SameLine(); ImGui::ProgressBar(progress);
ImGui::End();
}
/*
ImGui::SetNextWindowSize(ImVec2(1000, 400), ImGuiCond_FirstUseEver);
ImGui::SetNextWindowPos(ImVec2(680, 50), ImGuiCond_FirstUseEver);
{
ImGui::Begin("Ptolemy");
const int maxCount = 20;
static char reactionChar[maxCount][40];
static char gsSpinChar[maxCount][10];
static char KEAChar[maxCount][10];
static char orbitalChar[maxCount][10];
static char ExChar[maxCount][10];
static char jpiChar[maxCount][10];
static char potInChar[maxCount][10];
static char potOutChar[maxCount][10];
if( finsihed ){
if( ImPlot::BeginPlot("Plot", ImVec2(ImGui::GetWindowWidth()-20, 400)) ){
ImPlot::SetupLegend(ImPlotLocation_East, ImPlotLegendFlags_Outside);
ImPlot::SetupAxes("A"," Energy [MeV]");
ImPlot::SetupAxesLimits(ARange[0], ARange[1], -50, 0);
// ImPlot::SetupAxisLimitsConstraints(ImAxis_X1, ARange[0], ARange[1]);
// ImPlot::SetupAxisLimitsConstraints(ImAxis_Y1, floor(V0)*2, 0);
const double * kaka = AList.data();
for( int i = 0; i < orbitalStr.size(); i++){
vector<double> jaja;
for( int j = 0; j < AList.size(); j++) jaja.push_back(enLevelList[j][i]);
const double * haha = jaja.data();
ImPlot::PlotLine(orbitalStr[i].c_str(), kaka, haha, AList.size());
}
static int count = 1;
ImGui::InputInt("Number of states", &count, 1);
ImPlot::EndPlot();
const int nCol = 8;
if( ImGui::BeginTable("Input", nCol, ImGuiTableFlags_None) ){
{//Header
ImGui::TableNextColumn(); ImGui::Text("Reaction");
ImGui::TableNextColumn(); ImGui::Text("G.S. Spin");
ImGui::TableNextColumn(); ImGui::Text("Energy [MeV/u]");
ImGui::TableNextColumn(); ImGui::Text("Orbital");
ImGui::TableNextColumn(); ImGui::Text("Ex [MeV]");
ImGui::TableNextColumn(); ImGui::Text("J-pi");
ImGui::TableNextColumn(); ImGui::Text("Pot-in");
ImGui::TableNextColumn(); ImGui::Text("Pot-out");
ImGui::TableNextRow();
}
static bool saveDataClicked = false;
if( ImGui::Button("Download data") ){
// FILE * file_out;
// file_out = fopen("", "w+");
saveDataClicked = true;
}
if( saveDataClicked ) {
ImGui::Text("to be impletmented....");
}
for( int i = 0; i < count; i++){
ImGui::TableNextColumn(); ImGui::InputText(" ", reactionChar[i], 40) ;
ImGui::TableNextColumn(); ImGui::InputText(" ", gsSpinChar[i], 10) ;
ImGui::TableNextColumn(); ImGui::InputText(" ", KEAChar[i], 10) ;
ImGui::TableNextColumn(); ImGui::InputText(" ", orbitalChar[i], 10) ;
ImGui::TableNextColumn(); ImGui::InputText(" ", ExChar[i], 10) ;
ImGui::TableNextColumn(); ImGui::InputText(" ", jpiChar[i], 10) ;
ImGui::TableNextColumn(); ImGui::InputText(" ", potInChar[i], 10) ;
ImGui::TableNextColumn(); ImGui::InputText(" ", potOutChar[i], 10) ;
ImGui::TableNextRow();
}
ImGui::EndTable();
}
if( ImGui::Button("Calculate") ){
///==== create infile
///==== run Ptolemy
///==== extract d.s.c. from out file
///==== plot
}
ImGui::End();
}
*/
ImGui::Render();
int display_w, display_h;