Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
schroedinger_hydrogen.C
Go to the documentation of this file.
1/// \file
2/// \ingroup tutorial_graphics
3/// \notebook
4/// Plot the Amplitude of a Hydrogen Atom.
5///
6/// Visualize the Amplitude of a Hydrogen Atom in the n = 2, l = 0, m = 0 state.
7/// Demonstrates how TH2F can be used in Quantum Mechanics.
8///
9/// The formula for Hydrogen in this energy state is \f$ \psi_{200} = \frac{1}{4\sqrt{2\pi}a_0 ^{\frac{3}{2}}}(2-\frac{\sqrt{x^2+y^2}}{a_0})e^{-\frac{\sqrt{x^2+y^2}}{2a_0}} \f$
10///
11/// \macro_image
12/// \macro_code
13///
14/// \author Advait Dhingra
15
16#include <cmath>
17
18double WaveFunction(double x, double y) {
19 double r = sqrt(x *x + y*y);
20
21 double w = (1/pow((4*sqrt(2*TMath::Pi())* 1), 1.5)) * (2 - (r / 1)*pow(TMath::E(), (-1 * r)/2)); // Wavefunction formula for psi 2,0,0
22
23 return w*w; // Amplitude
24
25}
26
27void schroedinger_hydrogen() {
28 TH2F *h2D = new TH2F("Hydrogen Atom",
29 "Hydrogen in n = 2, l = 0, m = 0 state; Position in x direction; Position in y direction",
30 200, -10, 10, 200, -10, 10);
31
32 for (float i = -10; i < 10; i += 0.01) {
33 for (float j = -10; j < 10; j += 0.01) {
34 h2D->Fill(i, j, WaveFunction(i, j));
35 }
36 }
37
40
41 TCanvas *c1 = new TCanvas("c1", "Schroedinger's Hydrogen Atom", 750, 1500);
42 c1->Divide(1, 2);
43
44 auto c1_1 = c1->cd(1);
45 c1_1->SetRightMargin(0.14);
46 h2D->GetXaxis()->SetLabelSize(0.03);
47 h2D->GetYaxis()->SetLabelSize(0.03);
48 h2D->GetZaxis()->SetLabelSize(0.03);
49 h2D->SetContour(50);
50 h2D->Draw("colz");
51
52 TLatex *l = new TLatex(-10, -12.43, "The Electron is more likely to be found in the yellow areas and less likely to be found in the blue areas.");
53 l->SetTextFont(42);
54 l->SetTextSize(0.02);
55 l->Draw();
56
57 auto c1_2 = c1->cd(2);
58 c1_2->SetTheta(42.);
59
60 TH2D *h2Dc = (TH2D*)h2D->Clone();
61 h2Dc->SetTitle("3D view of probability amplitude;;");
62 h2Dc->Draw("surf2");
63}
@ kCividis
Definition TColor.h:128
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t r
R__EXTERN TStyle * gStyle
Definition TStyle.h:433
virtual void SetLabelSize(Float_t size=0.04)
Set size of axis labels.
Definition TAttAxis.cxx:203
The Canvas class.
Definition TCanvas.h:23
TAxis * GetZaxis()
Definition TH1.h:324
void SetTitle(const char *title) override
Change/set the title.
Definition TH1.cxx:6707
TAxis * GetXaxis()
Definition TH1.h:322
TAxis * GetYaxis()
Definition TH1.h:323
void Draw(Option_t *option="") override
Draw this histogram with options.
Definition TH1.cxx:3067
virtual void SetContour(Int_t nlevels, const Double_t *levels=nullptr)
Set the number and values of contour levels.
Definition TH1.cxx:8400
TObject * Clone(const char *newname="") const override
Make a complete copy of the underlying object.
Definition TH1.cxx:2734
2-D histogram with a double per channel (see TH1 documentation)}
Definition TH2.h:301
2-D histogram with a float per channel (see TH1 documentation)}
Definition TH2.h:258
Int_t Fill(Double_t) override
Invalid Fill method.
Definition TH2.cxx:347
To draw Mathematical Formula.
Definition TLatex.h:18
virtual void Draw(Option_t *option="")
Default Draw method for all objects.
Definition TObject.cxx:274
void SetOptStat(Int_t stat=1)
The type of information printed in the histogram statistics box can be selected via the parameter mod...
Definition TStyle.cxx:1636
void SetPalette(Int_t ncolors=kBird, Int_t *colors=nullptr, Float_t alpha=1.)
See TColor::SetPalette.
Definition TStyle.cxx:1884
RVec< PromoteTypes< T0, T1 > > pow(const T0 &x, const RVec< T1 > &v)
Definition RVec.hxx:1809
Double_t y[n]
Definition legend1.C:17
return c1
Definition legend1.C:41
Double_t x[n]
Definition legend1.C:17
VecExpr< UnaryOp< Sqrt< T >, VecExpr< A, T, D >, T >, T, D > sqrt(const VecExpr< A, T, D > &rhs)
constexpr Double_t E()
Base of natural log: .
Definition TMath.h:93
constexpr Double_t Pi()
Definition TMath.h:37
TLine l
Definition textangle.C:4