The ignorance of West Bengal Govt. about fundamental science

A proposal of a billion dollar science experimental facility to build a electron cyclotron in India was put forward in 2011 by the international committee headed by SINP (saha institute of nuclear physics,kolkata) a premier institute in this field.The machine Indian Synchrotron for Materials and Energy Research (ISMER) would be one of the largest in the world after USA,Germany,Japan,France. But only iff west bengal govt offered a 100 Acre land within the span of 50 kms from city and airport for easy transport for physicists across borders . But this is too much to ask for the government who believes in immediate benefit from investment and insists that it will be better to give the land for building AIIMS( medical school) than GOD-KNOWS-WHAT electron synchrotron. This chaps have not even heard of one before, although,  one similar technology is running right in front of their nose at VECC (variable energy cyclotron centre).

Bengal has been the hotspot for fundamental Physics education for long time. It has developed into one of the leading center for making of good scientists starting from Prof. Satyendra Nath Bose, the man who realized the famous “Bose particle”, a now fundamental part of physics to recently Dr. Ashok Sen for his outstanding research in fundamental string theory. What these two men have in common is they hailed from a place called the “intellectual” capital of India:West Bengal and have devoted their lives trying to understand the fundamental constituents of nature.

But the same Bengal, with due respect, is the home of another kind of intellectuals, the “wanna be” intellectuals, who quickly climb up ladder of high politics, who I guess cannot fathom the scale of development and progress a science project can build.

Coming to the points why this should be build:

1. Every nation must strive to achieve something hard/majestic, that is what sets them apart, remember them, just like Maharajahs are remembered by their wonderful monuments which sometimes took decades to complete, the modern day monuments are the Science and Technological marvels. Constructing this will take 10 years almost and employ thousands of people and also advance technology manifold.

2.Kolkata is home of many great scientist, and having this in homeland will make us even more proud.

3.Wiki gives the following uses for Science and Technology apart from fundamental physics:

Bengal have now been producing some software technicians like one in each hundred, only to brain drain to other state or even country doing some poor quality work for foreign demands. This project will employ enthusiast engineers in thousands. I have seen people from Bengal holding high posts in other states, intelligent people, scientists, artists fleeing this land, because of lack of opportunity, lack of cooperation from government, lack of understanding amongst those who take the decisions. Lets bring in a change, a change so that many people can really build something, a thing that will be remembered in history. Lets aim for long term benefit.

http://articles.timesofindia.indiatimes.com/2013-06-17/kolkata/40026493_1_sinp-deutsches-elektronen-synchrotron-saha-institute

http://www.telegraphindia.com/1130301/jsp/nation/story_16619819.jsp#.UgITxOGOWBs

http://www.telegraphindia.com/1130617/jsp/calcutta/story_17016076.jsp#.UgITzeGOWBs

2D Fast Poisson Solver

2D Fast Poisson Solver:

Introduction:

Here we will construct poisson solver for a 2D square boundary using a variation
fast fourier transform , dicrete sine transform.There are 4 types of discrete sine transform,you can read more about them in
wikipedia http://en.wikipedia.org/wiki/Discrete_sine_transform

We will solve the 2D electrostatic equation for square boundary,in

\epsilon_0=1 units

The poisson equation solver in 2D basically uses the matrix representation of the discrete poisson equation.
The poisson equation \nabla^2 u= \rho can be solved using fourier transform technique.The algorithm is basically
from
http://www.math.nus.edu.sg/~bao/teach/ma5233/lect10.ppt

This uses dst function(Type -I) from either fft.fftpack or a seperate
dst package I downloaded from internet called transforms

CODE:

In [2]:
from __future__ import division
import matplotlib.pyplot as plt
import numpy as np
from numpy.fft import fft,fft2,ifft2,ifft,irfft2,rfft2
import random as random
from mpl_toolkits.mplot3d import Axes3D

Define a function for discrete sin transform using fft algorithm,by symmetrizing the function x

In [3]:
def dst(x,axis=-1):
    """Discrete Sine Transform (DST-I)

    Implemented using 2(N+1)-point FFT
    xsym = r_[0,x,0,-x[::-1]]
    DST = (-imag(fft(xsym))/2)[1:(N+1)]

    adjusted to work over an arbitrary axis for entire n-dim array
    """
    n = len(x.shape)
    N = x.shape[axis]
    slices = [None]*3
    for k in range(3):
        slices[k] = []
        for j in range(n):
            slices[k].append(slice(None))
    newshape = list(x.shape)
    newshape[axis] = 2*(N+1)
    xsym = np.zeros(newshape,np.float)
    slices[0][axis] = slice(1,N+1)
    slices[1][axis] = slice(N+2,None)
    slices[2][axis] = slice(None,None,-1)
    for k in range(3):
        slices[k] = tuple(slices[k])
    xsym[slices[0]] = x
    xsym[slices[1]] = -x[slices[2]]
    DST = fft(xsym,axis=axis)
    #print xtilde
    return (-(DST.imag)/2)[slices[0]]

Define 2 dimensional DST, the idst is same as DST because of symmetry

In [4]:
def dst2(x,axes=(-1,-2)):
    return dst(dst(x,axis=axes[0]),axis=axes[1])

def idst2(x,axes=(-1,-2)):
    return dst(dst(x,axis=axes[0]),axis=axes[1])

Now define the poisson solver in 2D square boundary using the algorithm described above,this method can be
generalised to arbitary retangular boundary

In [5]:
def fft_poisson(f,h):

	m,n=f.shape
	f_bar=np.zeros([n,n])
	u_bar = f_bar			# make of the same shape
	u = u_bar

	f_bar=idst2(f)			# f_bar= fourier transform of f

	f_bar = f_bar * (2/n+1)**2  #Normalize
	#u_bar =np.zeros([n,n])
	pi=np.pi
	lam = np.arange(1,n+1)
	lam = -4/h**2 * (np.sin((lam*pi) / (2*(n + 1))))**2 #$compute $\lambda_x$
	#for rectangular domain add $lambda_y$
	for i in xrange(0,n):
		for j in xrange(0,n):
			u_bar[i,j] = (f_bar[i,j]) / (lam[i] + lam[j])

	u=dst2(u_bar)				#sine transform back
	u= u * (2/(n+1))**2			#normalize ,change for rectangular domain
	return u

Fix the initial parameters,domain of the equation etc.

In [6]:
#set bounds a,b,parameters
a = 0; b = 1;	
alpha=10				#alpha is grid points=2^alpha
n=2**alpha
L=b-a					#length of system

xe =np.linspace(a,b,n); 
ye =np.linspace(a,b,n);
x, y = np.meshgrid(xe,ye)
In [7]:
h=L/(n);			#size 
h2 = h**2;			#h squared
hx2=h2; 			#started with a cube,hx2=hy2
hy2=h2;
f=np.zeros([n,n]);

We put charges in the field.

Putting a random negetive charge in the field alog with a posituve charge

In [8]:
#initial conditions

#f[round(n/2):n,round(n/2):n]=1; #Initial condition
f[round((n+1)/2),round((n+1)/2-10)]=20
#f[round((n+1)/2),round((n+1)/2+10)]=-20
f[random.randint(0,n),random.randint(0,n)]=-10	#put a random charge
In [9]:
nx,ny=np.array(f.shape)-1		 #last index
U=np.zeros([n,n])

set boundary conditions

In [10]:
# BOUNDARY CONDITIONS
#set U(x,y)=g(x,y)
U[0,:]=0
U[nx,:]=0#5e-5
U[:,0]=0
U[:,ny]=0
In [11]:
##homogenize boundary condition

f[1,:]=f[1,:]+U[0,:]/hx2;
f[nx-1,:]=f[nx-1,:]+U[n-1,:]/hx2;

f[:,1]=f[:,1]+U[:,0]/hy2;
f[:,ny-1]=f[:,ny-1]+U[:,n-1]/hy2;

Solve the equation:

In [12]:
U=fft_poisson(f,h)
In [13]:
plt.figure()
#plt.imshow((U),cmap='hot')
plt.contour(U,50)
#plt.contour(f)

plt.show()
In [14]:
display(gcf())

In [15]:
imshow(U)
Out[15]:
<matplotlib.image.AxesImage at 0xb2b8a8c>
In [16]:
display(gcf())

 

In [ ]:

//

Hall Effect- a different outlook for seeing Hall coefficient

Hall effect is an extraordinary phenomena which involves generation of a subsidiary Electric field V_H when a electric field and magnetic field is applied to a conductor. You can get a basic introduction from wiki. If you are not familiar with it

But familiar Hall effect introduces it through a bound rectangular boundary, which limits the current destity in y direction to zero in steady state. In brief the Hall coefficient is defined as R_H=\frac{E_y}{J_x B}=\frac{1}{ne}

where ”J_x” is the current density of the carrier electrons, and E_y is the induced electric field,” n” is the carrier charge density and e is the carrier charge.

But consider the Hall effect in a 2D ring.The electric field applied in the ends if inner and outer ring and the magnetic field in the transverse direction.Then solving the Lorentz force equation along with Drude theory model (eqn 1)we get a 2 coupled equation,which can be solved directly

\frac{dp}{dt}=-e(E+v\times B)-\frac{p}{\tau} —–(1)

which gives from drude theory (see Ascroft Mermin)

\sigma_0 E_x=\omega_c \tau J_y+J_x ----(2) \\    \sigma_0 E_y=-\omega_c \tau J_x+J_y----(3)

Solving the two equations by multiplying the equations successively by p_x,p_y and adding/sutracting we get

\sigma_0 (E\times J)_z=\omega_c \tau J.J—-(4)

and

\sigma_0 (E. J)_{xy}= J.J—-(5)

where \sigma_0 (E\times J)_z=E_xJ_y-E_yJ_x—-(6)

\sigma_0 (E. J)_{xy}=J_xE_x+J_yE_y—-(7)

So we can define tan (\theta_H)=\omega_c\tau

which in common terms is the hall angle

where \omega_c=\frac{eH}{m}

Therefore the relation can be written in a more generally

\frac{tan (\theta_H)}{H}=\frac{\sigma}{ne} \\using \\ \sigma=\frac{ne^2\tau}{m}

The hall coefficient would then be defined by

R_H=\frac{tan (\theta_H)}{\sigma H}=\frac{1}{ne}

The angle between the net electric field and the net current is given by the angle tan (\theta)

This to me is a more general statement for Hall coefficient which reduces to the normal hall coefficient if we consider appropriate boundaries and noticing tan(theta) must be dimensionless and considering E_x<<E_y we get

R_H=\frac{E_y}{J_x H}=\frac{1}{ne}

which is the usual hall-coeficient.

Therefore the definition of Hall coefficient should be the tangent of Hall angle between the net Electric field and current in the conductor per unit applied magnetic field per unit conductivity.

Also this would help to resolve some issues regarding energy conservation in hall effect.

Suggestions are welcome.

Aritra Kundu

PS:here is an interesting article regarding hall effect in circular geometry http://www.jstor.org/stable/983911?seq=4

DL_POLY a quick user guide-Installation to first simulation in Molecular Dynamics

This summer I spend my time in Ecole Centrale Paris, trying to run a computer simulation on NPLIN.I will post about that later, but this was my first introduction to molecular dynamics simulation,  and there I faced many problems.So I thought to write up a bit, in case it might be of help for someone in future.I am posting this , but you sure can find help from the DISCO forums.And yes, read through the manual!!

In this, I explain how we can run DL_POLY in parallel.It really did a 30 minuite simulation in 30 seconds in my 4X2 core processor in laptop

Installation:
1.Download DL_POLY_4 from ftp://ftp.dl.ac.uk/ccp5/DL_POLY/
2.Obtain Liscense by mailing the authors for home/personal/educational use
3.Unzip,Untar,Unencrypt
4.In Linux install open mpi by typing
$ sudo apt-get install libopenmpi-dev openmpi-bin openmpi-doc
$ sudo apt-get install ssh
Generate keys and add to ssh
$ ssh-keygen -t dsa
$ cat id_dsa.pub >> authorized_keys

Test program:
#include <stdio.h>
#include <mpi.h>

int main(int argc, char** argv)

{
int myrank, nprocs;

MPI_Init(&argc, &argv);
MPI_Comm_size(MPI_COMM_WORLD, &nprocs);
MPI_Comm_rank(MPI_COMM_WORLD, &myrank);
for (;;){
printf(“Hello from processor %d of %d\n”, myrank, nprocs);
}
MPI_Finalize();
return 0;
}

Testing:
$ mpicc MPI_Hello.c -o MPI_Hello
$ mpiexec -n 5 MPI_Hello         #specify to run on 5 processors,press Ctrl+Z to exit
You can install htop to see if the process is really parallel [sudo apt-get install htop]
Compiling DL_POLY for Laptop Dell XPS 15z
Requires fortran to be installed.[sudo apt-get install gfortran]
1.In the DL_POLY directory>Build>copy “Makefile_MPI” to source directory and rename to “make”
2.In terminal type “sudo make hpc”
3.The process will create a DLPOLY.Z in the execute directory. Which also has a JAVA gui which can be used.(Personally I find using GUI a bit cumbersome and problematic,so I run using terminal)
4. Terminal >browse to DL_POLY execcute directory> type $./gui to invoke gui

Running sample scripts:
Some example and basic examples can be found online in
http://www.ccp5.ac.uk/DL_POLY/TUTORIAL/EXERCISES/

Although they used GUI, I am going to use a mixture of both

DLPOLY requires some necessary files to run a simulation,The CONTROL,CONFIG,FIELD files.I am going to explain them as we go,but first lets see if its working!
Download the files for EXERCISE 1 which finds out the phase transition for KCl.
The CONFIG file lists the system in brief,it gives the coordinates and velocities of all atoms in the system in initial condition.
The CONTROL file is important,it determins the thermodynamics of the system.The temperature, pressure, cutoff, and configuration like NVT,NPT etc.
FIELD defines any other field which cannot be characterised by concrete mathematical formulas

The CONFIG file uses its own self-consistent set of units.

physical quantity symbol unit value
time t0 1×10−12 seconds (picoseconds)
length l0 1×10−10 metres (Angstroms)
mass m0 1.66054×10−27 kilograms (amu)
charge q0 1.60218×10−19 Coulombs (electron charge)
energy E0=m0(l0/t0)2 1.66054×10−23 Joules = 10 J mol−1
pressure p0=E0/l30 1.66054×107 Pascal = 166.054 bar
Planck constant 6.35078E0t0
Boltzmann constant kB 0.83145E0/K

The others can be found from the user manual(Page 7 of DL_POLY_4_USER_MANUAL)

will learn the process through a process through the first example of DL_POLY guide.i.e. phase transition of KCl

http://www.ccp5.ac.uk/DL_POLY/TUTORIAL/EXERCISES/expmt1.html

This exercise studies a well-known phase transition in potassium chloride using constant pressure molecular dynamics. The objective is to develop the best practice in using such algorithms and to learn how phase transitions can be induced, detected and monitored in a simulation.

The CONFIG file:

DL_POLY: Potassium Chloride Test Case                                           #heading
2         3      2000    0.5000000000E-02

!levcfg(0-2,2 :position,velocity,forces included in file)

! imcon(periodic boundary condition:0,1,2,3,6)

!megatm(total number of particles)

!crystallographic information for periodic boundary(x,y,z component of a,b,c cell vector)

18.6960000000        0.0000000000        0.0000000000
0.0000000000       18.6960000000        0.0000000000
0.0000000000        0.0000000000       18.6960000000
!since levcfg = 2 each record will contain 3 data position,veocity,force
K+               1                            #atom name, index
-7.608595309        -7.897790000        -7.892053559    #position
1.056610291        -1.218664448         3.345828610    #velocity
-1979.558687         739.7961625         1027.996603    #force
K+               2
-4.678562643        -4.509858530        -7.813252872
-0.7078383910         1.220785303       -0.9440783885
212.7300659        -4290.721728         1363.866117
K+               3
-7.821414265        -4.635443539        -4.732164540
2.637614561        0.5778767520E-01    -1.704765568
150.3309776        -812.6932914         1429.413120
K+               4
-4.577735234        -7.808222424        -4.540282303
4.950411107         3.576253827        0.2276058039
-782.7309730         1322.553499        -205.1240756
Cl-              5
-7.810782981        -4.845250366        -7.815286768
1.643324971       -0.6479194620         3.393130560
407.2394357         1675.508460         748.2506938

The details about preparing this file with other packages(Athen) will be discussed
Next we prepare the CONTROL file:

**–begins here–<We will at first discuss only relevant parameters>–**
DL_POLY: Potassium Chloride Test Case    !heading
temperature      600.0                !system temperature in Kelvin
pressure         144.0                !system pressure
steps             4000                !maximum no of steps,2000 is enough
equilibration      500                  !time for the system to equlibrate
timestep         0.005
multiple             1
print               10
stats               10
scale               10
restart scale
rdf                 10             !write radial distribution function at 10 timescale interval
print rdf
cutoff             7.0                !define cutoff for interactions to 7 Angstorms
rvdw               7.0
delr               1.0
ewald precision  1.e-6
ensemble nst hoover 0.1 1.0            !Ensemble
trajectory 20 30 0                !write history file
job time          3600                !set job time
close time         100
finish                        !control file ends

FIELD FILE:

DL_POLY: Potassium Chloride Test Case

UNITS eV    !energy units

MOLECULES 1    ! number of molecule type

Potassium Chloride    !Molecule name

nummols 27        !number of molecules

atoms 8        ! no of atoms in this molecular type

K+           39.0983       0.994    4

Cl-          35.4530      -0.994    4

finish

VDW 3            !number of van der wall interactions

K+      K+      buck      3796.9      0.2603      124.90

Cl-     Cl-     buck      1227.2      0.3214      124.90

K+      Cl-     buck      4117.9      0.3048      0.0000

CLOSE

RUNNING DLPOLY:

To run the program make sure the CONFIG,CONTROL,FIELD files are there in the execute folder, from terminal invoke dlpoly

./dlpoly.z

 

This will run the simulation and write the output files

ANALYZING OUTPUTS:

 

  1. We will visualize the resulting system.This can be done through GUI by plotting the REVCON file(which is the system at the end of simulation)
    You can also visualize the CONFIG file.
    Another better alternative of visualizing the files is to convert the CONFIG file to XYZ files and plotting them on VMD

  2. Checking the Diffusion coefficient in OUTPUT file, if it is less than 10-5then it is liquid state.A value of zero suggests it to be solid state.

  3. Now we will analyze the output from the simulation using the RDFDAT file that was created.From GUI plot the RDF -> select the atom types K+ and Cl-

    This will extract the RDF for K+ and Cl- from the RDFDAT file and create a separate XY file.For various simulation parameters you can create different RDF files and merge them together and plot so that you can compare where the system is heading to and draw good conclusions.I have written a small script in python to do that.

# Import modules

import sys

import matplotlib.pyplot as plt

import numpy as np

”’

to be invoked as

plot_rdf <option={-m,-p}> <no of files> <list of files>

option -m merges the given number of files using single X cordinates and

y coordinate.Option -p plots the given file

eg

plot_rdf -m 2 RDF_100.XY RDF_14.XY

merges the two files

”’

#import pylab as pl

option=sys.argv[1]

no_of_options=int(sys.argv[2])

index=3

l_index=index+no_of_options

#process data

if(option==’-m’):

   data=np.loadtxt(sys.argv[index])

   shape_data=data.shape

   while(index<l_index-1):

       index=index+1

       data2=np.loadtxt(sys.argv[index])[:,1].reshape(shape_data[0],1)

       data=np.append(data,data2,1)

   f=open(‘RDF_Merged.XY’,’w’)

#TODO: write names of field and parameters in first column

   for row in data:

       X=”

       for element in row:

           X=X+’ ‘+str(element)

       f.write(X+’\n’)

   f.close()

elif (option==’-p’):

   # Open the data

   datafile = sys.argv[index]

   #f = open(datafile, ‘r’)

   data = np.loadtxt(datafile)

   shape_data=data.shape

   i=1

   #plot all Y’s

   while(i<shape_data[1]):

       plt.plot(data[:,0],data[:,i])

       i=i+1

   # Get axis labels

   xaxis_label = raw_input(‘X-axis label > ‘)

   yaxis_label = raw_input(‘Y-axis label > ‘)

   # Show the plot

   plt.ylabel(yaxis_label)

   plt.xlabel(xaxis_label)

   plt.show()

#Stop

The RDF is the measure of how the probability deviates from an uniform distribution.From the radial distribution function a number of parameters can be determined.We can also identify the phase transition from the change in pattern in the RDF

Picture stolen from Jacek Dziedzic,Gdansk University of Technology,Computational nanotechnology

google ngram observation

I was experimenting with the Google Ngram
Its a wonderful tool that i came across after watching the video on “what we learned from 5 million books”from TedX

For a time now i had wondered when this new age saga for Aliens begin,and the term UFO and Flying saucer are closely related to the controversy of aliens ,So I did a quick Ngram search on UFO and aliens.The results pointed to interesting results.The use of these words started round about 1950’s,but critically the word UFO was more convinient and ahd the mysterious feeing which made it take over the word flying saucer.The heights of the UFO alien visit was in 2000’s when (both the words show a peak) .Well I still guess aliens are out there,
Daleks will be coming to earth one day,and Dr.Who where are you????

SETI institute should be funded for better research for time travel and communication with ET’s

Toward Politics:
The two main groups and thought of interest in this world and its shaping to its current state have been the capitalist and the communist.The Ngram tells us communism was old concept than the cpitalist at least in words dating back at least to 1880.it hasd a slow rise among mass till in 1960’s  it achieved a peak, owing to the fast changing culture and socity in russia china and india
http://www.youtube.com/watch?v=D-ExLu2ilKE
Where as The capitaisit community was at peak in US round about 1940’s  during the war ,but with the rise of communism in 60’s the interest shifts to communusm.
But soon america became overpowerful and decied upon capitalism to bring forth power and trade .
The current trend seems a fast declining capitalism compared to communism.
Will the age of communism rise again?
Viva La revolution!

Ah ha here is another interesting stuff,cinema and movie..which came first and what is more common ?
surprisingly both were born round about the birth of film industry round about 1910.
The number of movies produced after 1960’s became large and also by the invention of color telivision the media was moved out to masses.
and well movie won with cinema..which was more awkward term more reserved for cinematographer. while movie was a sweet,short ,birsk word ready to be used in a date and with friends…
movies on rise..with time and pride…

Communism on its heights roun about 60’s!!


History of fuck:
In its modern sense it started it social outrage from little before 1960’sand then on rise ever…

An embarrassing comparison for google..yahoo seemed more written about till 2008..yes of course it should have been..any doubt????

Graviton research:
80’s!!!

Billion vs trillion:
I guess I was thinking in terms of money.we see with time more is being taked about heigher numbers.Are we heading towards society of more rich and profit?

Tracking history :
Discovery of pluto. round about 1940+-5 years for research to take place.Wiki says 1930’s wow close!!!