sample.orbital.density.mcmc.parallel: MCMC sample orbital density in separate processes

In npetraco/che302r: Handy functions for CHEM 302, Physical Chemistry II

MCMC sample orbital density in separate processes

Description

Parallel version of sample.orbital.density.mcmc

Usage

sample.orbital.density.mcmc.parallel(
  nqn,
  lqn,
  mqn,
  orb.func.info,
  printQ = FALSE,
  num.processes = 1
)

Arguments

nqn	n quantum number
lqn	l quantum number
mqn	m quantum number
orb.func.info	Orbital pdfs info generated by splined.sqorb.functs()

Details

Sample the wave function's density and plot corresponding r,theta, phi (x,y,z) values. Uses MCMC. Do each sampling chain in parallel. Under the hood, this is a VERY stripped down wrapper (to a wrapper: sample.density3b) to run MCMC. CHECK THAT YOU GET ALL THE NODES IN YOUR WAVE FUNCTIONS! MCMC by itself can easily miss them.

Value

orbital sample matrix (num.sample-rows by 3-columns)

References

XXXXX

Examples

library(che302r)
library(parallel)
library(doSNOW)
library(rgl)

#Pick an orbital:
n <- 1
l <- 0
m <- (0)

# Create squared orbital functions for sampling
orb.sq.funcs  <- splined.sqorb.functs(n, l, m, r.max = 320, num.knots = 1000)
rax           <- orb.sq.funcs$r.axis
R.sq          <- orb.sq.funcs$R.sq.func
plot(rax, rax^2 * R.sq(rax), typ="l")

#Sample the orbital density:
detectCores() # Number of processors available

#density.samples <- sample.orbital.density.ais(n, l, m, printQ = T, orb.func.info = orb.sq.funcs)
#density.samples <- sample.orbital.density.mcmc(n, l, m, printQ = T, orb.func.info = orb.sq.funcs)
#density.samples <- sample.orbital.density.ais.parallel(n, l, m, printQ = T, orb.func.info = orb.sq.funcs, num.processes = 2)
density.samples <- sample.orbital.density.mcmc.parallel(n, l, m, printQ = T, orb.func.info = orb.sq.funcs, num.processes = 4)

# Sample in cartesian coordinates
x               <- density.samples[,1]
y               <- density.samples[,2]
z               <- density.samples[,3]

#Plot the sampled orbital density:
#2D
plot(x,y)
plot(x,z)
plot(y,z)

#3D
plot3d(x,y,z,type="s",radius=1.5,xlab="x",ylab="y",zlab="z")

# 3D slices
plot3d(x,y,rep(0,length(z)),type="s",radius=0.5,xlab="x",ylab="y",zlab="z")
plot3d(x,rep(0,length(y)),z,type="s",radius=0.5,xlab="x",ylab="y",zlab="z")
plot3d(rep(0,length(x)), y, z,type="s",radius=0.5,xlab="x",ylab="y",zlab="z")

npetraco/che302r documentation built on April 17, 2025, 10:34 p.m.