Future Investment Trends in India: NIFTY 50 and Sectoral Analysis 2024


The content contains detailed insights and analysis on the Indian market, particularly focusing on the NIFTY 50 index and various sectors like Financial, Auto, IT, Midcap, Smallcap, and SENSEX. The reflections and predictions provide an interesting perspective on potential investment opportunities and market trends.

This post is mainly a self reflection post so that I do not forget what I have done from Macro analysis in Indian market. If any reader is curious of the details, feel free to message me. This is not a financial advice.

The analysis is based on the following observation: The timescale it takes to double an initial investment in various markets are different. Here I will focus on the Indian market. Specifically the analysis that I will show here is very simple-minded, approximate and is intended to serve as a poor-mans investment guide. However with appropriate mathematical tools one can make more precise decisions (which I will not publish here). This is meant as a rough guide for the future, but not a financial advice.

The NIFTY 50 (^NSEI) is a well diversified 50 stock index and it represent important sectors of the economy. The timescale for doubling is approximately 6-7 years as can be seen from the graph: it crossed 6000 mark in 2014 while the 12000 mark on 2020. This correlates to the fact that the Private equity flows mainly from the US in the market doubles every 7 years. Although there are some reasons to argue that the markets should move away from PE flows to FII flows, but that will take some time to occur. In the current scenerio, the price between 10000 and 20000 has taken ~3.5 years. It is safe to assume that the NSEI index fund will hover around 15000 to 30000 in the next 3.5 years ,i.e. upto 2027. One could use swing trading to make money out of this market in this timeframe. This is inline with the predictions of FII market inflow and decadal median returns.

The finer analysis is which sectors are supposed to grow? As the NIFTY index constituents change over time, predicting possible constituents gives an advantage over broader market. By looking at the FII inflows, we see the largest inflows in this month are in Financial sector, oil gas sector, IT and Automobile. The smallest inflows are in textiles, paper and media industry.

Now let us do a timescale analysis of these sectors,

FINANCIAL SECTOR:

The Finance sector has a doubling timescale of around 5 years (2013-2017). In the current cycle of the movement from 10000 to 20000, the cycle started in 2018 and ended at ~70% in 2022 (5 year cycle). It is well due now, in the next 5 years, the index will reach 40000. The third figure confirms continuation of slow growth and no possible immediate catastrophe nearby. The price will move to around 37000 where it will face a strong rejection.

AUTO SECTOR: The doubling timescale is 4 years (2012-2016). The price is expected to grow to 30000 where it will face a short resistance before continuing to 40000 by 2027-28.

IT SECTOR: The timescale of doubling is 6 years (2012-2018)

This is an interesting case, the 10000-20000 doubling should have happened from 2014 in 2020, but after the mini-crash around 15000, the market grew faster to reach 20000 in 2021. If we have the same timescale, the transition from 20000 to 40000 should occur around 2026-2027. The prices are already expensive for 2024, hence, we expect the price to decrease in 2024-2025 and then again increase in 2026.

MIDCAP INDEX: Doubling time is 7 years. Which makes the current prices extremely high. 20000-40000 should occur in 2028, it is expected that the midcap index will go down to 20000 level again.

SMALLCAP: The intrinsic doubling timescale is again 7 years. Current valuations tell me that the market already crossed the doubling time in 3.5 years. We expect a drawdown of the small cap sector to 20000 levels in coming years.

To wrap up, lets study the SENSEX: The double timescale is a bit more from 7-8 years. In this case, the price is a bit overvalued as the index should be at 80000 around 2027.

Summary: To conclude, the Auto sector and the Financial sector looks promising in the next 4 years. The Oil and Gas sector will also perform well, but I could not find an index related to it. Choose from the constituents of the Auto index, Oil and Gas index, Financial index and invest in only the large cap stocks. Another index which is at lows is the Media index, however, it is a relatively new domain and our methods are not reliable. This is the era of rich getting richer. Cheers and enjoy the ride. Not a financial advice. Python notebook:https://colab.research.google.com/drive/17wiML0tWF25aIXBJv6K-S3H01dCCN20Y?usp=sharing

Financial Learnings

Hi, I am bringing back to life this blog after many years.

The idea here would be to document some of my thoughts about my learnings from probability and statistics and applying it to financial learnings.

Occasionally I will predict the global structure of the world and some ramblings of where the world is heading towards.

Getting to know strings better in Python!

Hey welcome back!
Solutions to previous post questions:
Exercise: 1. define a function which does the following
print alternate letters
print the leftout letters
example : if string = ABCDEF # output : ACE,BDF # Hints :string[0::2]

In [108]:
def alternate(string):
    print string[0::2] # print alternate letters starting from 0
    print string[1::2] #print alternate letters starting from 1
In [109]:
alternate("ABCD123")#call the above "function"
AC13
BD2

In [110]:
alternate("I am a good girl")
Ia  odgr
 mago il

Ex2: Writing in reverse

In [112]:
"Hello World"[::-1] #last index -1 reverses the order of the sentence! basically it reads,
#start from 0 to end by going negative steps
Out[112]:
'dlroW olleH'

So we have seen upto now how we can define strings and access its letters in them. Now things are going to get interesting. We are now going to use “functions” or commands that already the computer knows to manipulate a string and make it do things. But first we need to know what the computer understands right? It might be a bit boring, and I am not going to give it to you all at once, but through tits and bits and you can find help from google as always! So lets define a string and begin playing with it a bit. Lets see how many ways we can define a string in python

In [1]:
string = "I am a String"
character='A'
multiple_string = """I am multiple line
                    String !
                    I am line 2!
                    Hey I am line3
                  """
In [2]:
string
Out[2]:
'I am a String'
In [3]:
character
Out[3]:
'A'
In [4]:
multiple_string
Out[4]:
'I am multiple line\n                    String !\n                    I am line 2!\n                    Hey I am line3\n                  '

O wait, what are those \n? Well computer understands some more letters which humans have no “sound” for. Example, “\n” stands for an “newline”(n for new) , it is basically equivalent of pressing enter! Similarly, “\t” is a tabular space! So in the above example a multiple line string is defined in triple quotes and it automatically takes in “newline” where ever you press “enter”

Ok now lets define some operation on strings. You can “add” two strings just as numerical things by “+” sign.(beware “-” does not work!) you can even multiply a sentence by a number, and it adds the same copy of itself , lets see some example below

In [8]:
string+string
Out[8]:
'I am a StringI am a String'
In [43]:
print "tabs\t"*4
# (the hash is for comment, nothing after 
#this is seen as a command its for you )see the \t after 'tabs'
print "tabs"*4# see the difference. *4 adds
#the sentence 4 times!
tabs	tabs	tabs	tabs	
tabstabstabstabs

Well now strings are immutable, i.e. you cannot change a string that you have defined.So any new changes you make must be saved under a new name to start with, you can also “reassign” the “changed” string to the original! Examples

In [47]:
str1  = "String1"
str1[:-1]# all but the last letter(called slicing)
Out[47]:
'String'
In [48]:
str1 #str1 remains the same
Out[48]:
'String1'
In [52]:
str = str1[:-1]#save it to a new variable 
In [55]:
str
Out[55]:
'String'
In [56]:
str1
Out[56]:
'String1'
In [57]:
str1 = str1[:-1] # or reassign to same variable name ,
#when you dont need the previous sentence anymore in future
In [59]:
str1# str1 is changed to "String"
Out[59]:
'String'

Well now we will digress a little and do some numbers and then come back and use “strings” with a lot more power!
We have seen how we can access each letters of the sentence by the positions. How now about we manually do this scan through each letters and do something with it. Lets first do a simple example, how about counting the number of ‘e’ in a sentence. (BTW, ‘e’ is the most frequently used letter in English)

In [60]:
string  = "Count the number of eee's in me."
In [62]:
length = len(string) # number of "letters", including space and punctuations in the sentence!save it to variable length
In [63]:
length
Out[63]:
32
In [69]:
count_e = 0 # this is my counter, if i find 'e' then i will add 1 to this
for i in xrange(length): # this is called a "loop",  i is like your index finger running across the length of the sentence
    
    if(string[i]=='e'): # check 'if' the letter at ith index is e. "==" checks if it is "equal". single "=" assigns value
        count_e = count_e+1 #if yes add 1 to count_e variable


 #Look carefully how I gave a "tab" in front of if. This is necessary and very important in python loops, to tell what are the 
 # set of instructions one will follow iteratively. Similarly for "if" we need to add in another "tab" to make it clear what are
 # the statements one need to follow if the statement is true!
In [70]:
count_e
Out[70]:
6

Now we will try different functions that are there in python, to do various jobs.BTW caps are important. meaning variable name
A = 2 and a = 2 are different!

In [73]:
string
Out[73]:
"Count the number of eee's in me."
In [79]:
string.count('e') # do the same thing (count no of occurrences)
                #   as above but this time it is inbuilt
Out[79]:
6
In [86]:
string.find('e')# find the "first" occurance of the letter in the string
Out[86]:
8
In [88]:
string.lower()# convert everthing to lower case
Out[88]:
"count the number of eee's in me."
In [92]:
string.upper()# convert everthing to upper
Out[92]:
"COUNT THE NUMBER OF EEE'S IN ME."
In [91]:
string.replace('e','_')# replace all occurrences with some other values
Out[91]:
"Count th_ numb_r of ___'s in m_."
In [100]:
"   abcd  ".strip()#removes "spaces" from lead or end
Out[100]:
'abcd'
In [103]:
string.split(' ')#split sentence into words by "space" = " "
Out[103]:
['Count', 'the', 'number', 'of', "eee's", 'in', 'me.']

Now let us count the total number of vowels in a sentence

In [117]:
sentence = "the quick brown fox jumped over the lazy dog" #define a sentence
length = len(sentence) # length of sentence(no of characters)
count_vowels = 0 # counter for vowels

for i in xrange(length): # same old indexing, get yours hands out!
    
    letter  = sentence[i] #current letter in the sentence
    if(letter =='a' or letter=='e' or letter=='i' or letter=='o' or letter=='u'):
        #check if the letter is "a" or "e" so on....
        count_vowels = count_vowels+1

print "Number of vowels:\t",count_vowels
print "Number of consonants:\t",length - count_vowels
Number of vowels:	12
Number of consonants:	32

Exercise:
A palindrome is a sentence, which when reversed reads the same as its original, other than punctuations and spaces. e.g.
“Was it a car or a cat I saw”
“Rise to vote Sir”

palindromecomics

Write a function to test if a sentence is a palindrome.
If it is a palindrome print “it is a palindrome”
else print “it is not a palindrome”

tips:
def isPalindrome(string):
#remove spaces (hint : strip)
#find reverse String
#check for equality of strings
print “it is a palindrome”

#call the function
isPalindrome(“Rise to vote sir”)

In []:

Installing python and your first program

Time has come to teach your computer to speak “Python”

Firstly, the choice of Operating system:

So what do you use?

Windows, Linux?

I would personally prefer Linux, as it is open source and that is what I primarily use for my work. Although it really does not matter once you teach your computer the language, it can understand what you say as well on windows or linux.

But hold on! I will be committing a big mistake if I don’t refer to the huge support group and documentation available online. I will not be extensive over these set of tutorials, but there is a huge crowd of people who is ready to help you always at stackoverflow.com and the Python Docs

Next, I would prefer one to use the userfrienidly IPython, but is it it better to go to the shop once rather than going 1000 times, here comes the Enthought distribution , select the operating system and the 32/64 bit system( which you can find from system properties) and download the package ,and install it! (Just do next next next  :p for windows). For linux its a bit tricky, you will download a “***.sh” file. save it somewhere, say downloads. Open terminal (I am assuming you have Ubuntu) with Ctrl+Alt+T( remember CAT!) , hopefully it will open in Desktop, change directory by command

cd ../Downloads

chmod +x “***.sh”

sudo ./***.sh

[Give your password]

This should be starting the installation process, and would ask for where you want to install just type “Y” or do as it suggests 🙂

#————————————————————————————————————————————————#

Once you are done installing , its time to start up the software,

In windows open “command” (start -> run->type “cmd”), a black box should pop up.

In Linux do the Ctrl+Alt+T thing again, and the lil box should show up.

Ok so there are few things you need to “know” to make things happen, its a repetitive thing, these are the few steps:

ok now type in:

ipython notebook –pylab [for a white good looking interface,opens up in browser] or

ipython –pylab [for nerdy looking ones]

I will probably follow notebook here, because its easier to blog in that 🙂

Once it fires up there is an option (New Notebook)

notebook

click on that, and here you have the first terminal,

 

Here is my first code, well as is the rule lets begin by "HelloWorld", its a computer speaking to human through human interface. Shift+Enter is the basic go to next line(Evaluate this line). So everything in this notebook is evaluated on the go, its like an advanced calculator.

In [2]:
print "Hello World"
Hello World

Yay , my computer responded back with "Hello World". Lets test if it can sum my grocery bill

In [3]:
print 212+190+40+15+64
521

Oh good, you see how easy a calculator it is! I can just sum eveything at one go!

Ok now remeber I told you to break up simple tasks, well lets do the hello world with "functions" now,( this is a small digresson, if you find it complicated, skip it, you will get to know its advantages in future)

In [4]:
def HW(string):
    print string
In [5]:
HW("Hello World")
Hello World

This will require some explation so, every small task must have some "name", here "HW" is such a name, "def" stands for defination! so in natural language it is equivalent to saying "define a task HW, which takes in a sentence and prints it out!" And then call it!

Ok lets now begin with some basics of text. So firstly note that "text" and numbers are different! And computers mostly understand numbers. Firstly for everything there is a name! even for a sentence!

In [6]:
sen = "AbCd123"
In [7]:
print sen # call sen by its name and see what is there # is a comment
AbCd123

  now a string is a collection of letters, so if say you wanted to know the 3rd letter of the sentence you would do something like:
In [8]:
sen[2] # notice for 3rd character we write 2, this is because the numbering in this language starts from 0
Out[8]:
'C'

This is the concept of an "list" , which is an collection of objects, here that is letters. In general an list is defined with brackets, like

In [9]:
LList = ["A","b","C","d","1","2","3"]
In [12]:
print LList[2]," ",LList[0]
C   A

one can also "slice" the array, that is tell the computer to return all letters after a particular letter, the syntax for that is simple like

In [14]:
sen[2:] #return everything else from 2nd position
Out[14]:
'Cd123'
In [16]:
sen[:2]#return everything else till 2nd position
Out[16]:
'Ab'
In [18]:
sen[::2]#return everything with a gap of 2, i.e. 1st letter then 3rd letter then 5th and so on
Out[18]:
'AC13'

In general,

In [20]:
sen[2:-1:2]#start from pos 2 and till the last, and jump by 2 steps in between
Out[20]:
'C1'

Exercise: 1. define a function which does the following

def stringP(string): # print alternate letters # print the leftout letters # example : if string = ABCDEF # output : ACE,BDF # Hints :string[0::2]

call the function with your favourite coute

Exercise 2: define a code to write the string in reverse, that is for input of Sen, we get '321dCbA'

In []:

Learning to speak to computers!

Good to be back! Forgive me of my sins for I was occupied with interests of my own. 🙂

So apparently my girlfriend is now studying Linguistics and programming can become handy for large dataset analysis and automating some stuffs, so she (who never programmed anything ever) asked me to teach “programming”. I thought for a while and thought it would be fun to teach someone programming from scratch, think about it, it is a “language” too, albeit a “written” one, which is quite similar to the set of instructions that we are accustomed to process.

Someone writes to you,

Find the vowels in the sentence:   “Programming is simple conversation with the machine!”

The word “Find” is a trigger your brain processes to “narrow” down your “attention” to your prior knowledge of vowels, namely “AEIOU”. You can do the same thing to a computer, using its language! Now as humans have developed different language, you can talk to computer using various “languages” which it can understand.

BTW computers mother toungue is just True(1) and Flase(0).

I have been a great fan of the language “Python”, because of its simplicity and user friendliness, so I will stick to that. The main aim of these set of blog tutorials will be 2 fold. Recently I have been interested in NLP, which is Natural Language Processing, and its applications, these blog posts will apart from doing the basics of programming, move on to some basic language processing tools, which will be my own study journal. I might think to introduce some maths later on in this course for understanding some statistical properties of language.

So are you ready to talk to computer and explore the fascinating world?

Lets begin with the analogies, so that it isnt too “unfamiliar”. Your friend who speaks only German says to you :

Schreib “back”

Now I cant understand what Schreib means! because I dont know the language! Similarly a computer would not know any language, but you can teach him by installing a “software” ! I will write about Installing python in your computer in the next post, but till then hang on. So the basic thing to know speak to computer is to know its language. Specifically forget the adjectives there are only Names(nouns) and Commands(verbs) :p

BTW did you know ?

Guido van Rossum, the creator of the Python language, named the language after the BBC show “Monty Python’s Flying Circus”. He doesn’t particularly like snakes that kill animals for food by winding their long bodies around them and crushing them.

Now it is important that you learn a language by practicing it! So it applies same here, the next set of tutorials will be on my view about how one should see programming, well since i have said that, I believe that we are fast evolving into a world where humans and technology will be inseparable, and programming will be an intrinsic part of everyone’s life, there will be pet cyborgs and humanoids to make your life easier, so you better learn how to control them.

Ok, enough digression, so lets get back to our second example, Schreib in german means write, to do the same thing to a computer you can just say

print “HelloWorld”

and as you might have imagined, it does exactly as it is told 🙂 We will get to these details later. To have a formal understanding on how programming works, you must be really clear on what you want to do, and how exactly you want to do it!

This is the first and last thing one must worry about. Like if I said, “you Good  are ” (well its a odd fact that you could still make out the meaning, but believe me “computers” take things quite “Literally”), this doesnt make sense because  what exactly I want to ask is “Are you Good?” So structure is important in program too. Again in the first example, how would you count all vowels in a sentence?

1. recall what vowels are! {A,E,I,O,U} , wow !

2. Read the sentence “Programming is simple conversation with the machine!” letter by letter and note down vowels from it!

3.Naievely, you can also take a copy and write 1+1+1….. as soon as you encounter a vowel in the sentence.

4. To be sane enough, start reading from left to right!

Simple isnt it? You see how I broke a “Task” to smaller steps of what needs to be done. This is an important concept.

Now first time my teacher would teach me how to find vowels like this, but the next time she says me to find vowels, its a “Task” we should already know. This is the concept of a function. It is always a good idea to learn to define new short “verbs” for actions which you do often, say you work in a brick Kline, and the owner asks you “to collect mud, shape it and then bake it”  , would it be easier for him to “invent” a single word  for the whole process say “comisab” or something! and everyone just understands what the owner is saying.

So I will end this topic by summarizing 2 things what i think is important for writing a program:

1. Manually first do what you want to instruct to your computer in pen and paper. Its like a translation, first your need to know your mother language!

Did you Know? Programme comes from word Pro meaning before, and graphe meaning write, which probably sums up,what you need to do before writing your first code!Think about it before writing!

2.Try to break the process to as simple as possible and group them as small “tasks”.

Next Tutorial:

Installing python in Linux machine/Windows. Editors

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

Fast Fourier Transform (FFT) examples

The Discrete Fourier Transform(DFT) is defined by X_k = \sum_{n=0}^{N-1} x_n \cdot e^{-i 2 \pi \frac{k}{N} n} latex and the inverse fourier transform is defined as x_n = \frac{1}{N} \sum_{k=0}^{N-1} X_k \cdot e^{+i 2 \pi \frac{k}{N} n} latex

In [2]:
from __future__ import division
import numpy as np
import matplotlib.pyplot as plt
import time

'''
	Define two functions to compute fourier transforms using two different methods
	1.	dft(x,inverse=False,timer=False): for INVERSE: inverse=true
						  for timed: timer=true
	2.	timedfft(x,inverse=False,timer):for INVERSE: inverse=true
						  for timed: timer=true
	3.	func(x): defines function at x
'''

def func(x):
    f=10.0 
    val=np.sin(2*np.pi*f*x)
    return val

'''
   	Computes the discrete fourier transform /inverse fourier transform
	using basic summation convention.
	OPTIONS:
		inverse: computes inverse fourier transform
		timer:	 prints time required for computation
'''
def dft(x, inverse = False, timer = False) :

    t = time.clock()							#get current time to calculate the computation time
    N = len(x)
    if inverse:
        inv=1
    else:
        inv=-1

    X =[0] * N								#make an array of size N

    for k in xrange(N) :						#compute the required sum
        for n in xrange(N) :
            X[k] += x[n] * np.e**(inv * 2j * np.pi * k * n / N)
        if inverse :
            X[k] /= N

    t = time.clock() - t						#get delta time
    if timer :
        print "Computed DFT of size",N,"in",t,"sec."
    return X

'''
	fourier transform using built in fast fourier method
'''

def timedfft(x,inverse=False,timer=False):

	t=time.clock()
	if not inverse:X=np.fft.fft(x)
	else:	X=np.fft.ifft(x)
	t=time.clock()-t
	if timer: print "Computed DFT using fft of size",len(x),"in",t,"sec."
	return X

'''
	Computing part: 1. Define a function over a interval
			2. Compute fourier transform using the two methods
			3. Compare time
'''

#fix sample size, this gives the the rate of sampling of the function
sample_size=2**8
L=5.0
x=np.linspace(0,L,sample_size)

#there are four plots 1.original signal 2.fourier transformed signal 3.fast fourier transformed signal 4.Reconstructed signal

plt.subplot(411)
plt.xlabel('t')
plt.ylabel('f($t$)')
plt.xlim([0,L/10])
signal=func(x)
plt.plot(x,signal,label='original signal')

'''
				In [20]: run fourier_basic.py
				Computed DFT using fft of size 256 in 0.0 sec.
				Computed DFT of size 256 in 0.48 sec.
				Computed DFT of size 256 in 0.16 sec.

				In [21]: run fourier_basic.py
				Computed DFT using fft of size 4096 in 0.0 sec.
				Computed DFT of size 4096 in 123.94 sec.
				Computed DFT of size 4096 in 43.18 sec.
'''
#compute timed fft
fourierfft=timedfft(signal,timer=True)

#caluclate discrete fourier transform using normal method
fourier =dft(signal,timer=True)

#get the frequencies
n=x.size
d=x[1]-x[0]
freqsfft=np.fft.fftfreq(n,d)

#get positive index
indices=np.where(freqsfft>=0)

#plot
plt.subplot(412)
plt.xlabel('freq')
plt.ylabel('f($\omega$)')
plt.plot(freqsfft,fourierfft,label='FFT')

plt.subplot(413)
plt.xlabel('freq')
plt.ylabel('f($\omega$)')
plt.plot(freqsfft,fourier,color='green',label='Without FFT')

#inverse fourier transform using dft
invfourier=(dft(fourier,inverse=True,timer=True))

plt.subplot(414)
plt.xlabel('t')
plt.ylabel('f($t$)')
plt.xlim([0,L/10])
plt.plot(x,invfourier,'+-',color='red')
display(gcf())
Computed DFT using fft of size 256 in 0.0 sec.
Computed DFT of size 256 in 0.48 sec.
Computed DFT of size 256 in 0.16 sec.

/numpy/core/numeric.py:235: ComplexWarning: Casting complex values to real discards the imaginary part
  return array(a, dtype, copy=False, order=order)

N.B.: The fourier transform of the function is symmertric about 0.
Thus it suffices to display only the positive half,
as the negetive half gives us no more information

In [3]:
clf()
plt.subplot(111)
pi=np.pi
sampled=5*np.sin(2*pi*x)+np.sin(2*pi*10*x)
fourier=np.fft.fft(sampled)
#if you want to print only the positive frequencies
index=np.where(freqsfft>=0)
plt.plot(freqsfft[index],fourier[index])
display(gcf())

In [39]:
clf()
sampled=np.exp(-x**2)
fourier=np.fft.fft(sampled)
plt.plot(freqsfft,fourier)
display(gcf())

In [40]:
clf()
from scipy.special import jn
sampled=jn(2,x)
fourier=np.fft.fft(sampled)
plt.plot(freqsfft,fourier)
display(gcf())

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 [ ]:

//