Skip to content

Reimann Integration with Mayavi2 and traits

January 2, 2011

I have already blogged Reimann Integration using Mayavi2 here https://dhastha.wordpress.com/2010/12/26/reimann-integration-visualization-using-mayavi2/. I have upgraded this script with traits and Mayavi2. We can change radius of circle and number of slice to divide the circle as our wish, and we can see various result instantly on scene.

from enthought.traits.api import HasTraits, Range, Instance, on_trait_change
from enthought.traits.ui.api import View, Item, HGroup, HSplit, Group
from enthought.tvtk.pyface.scene_editor import SceneEditor
from enthought.mayavi.tools.mlab_scene_model import MlabSceneModel
from enthought.mayavi.core.ui.mayavi_scene import MayaviScene
from enthought.tvtk.pyface.scene import Scene
import numpy as np
from enthought.mayavi import mlab

x0,y0,z0=0,0,0

def circlepoints(r, n):
    a=n*100
    t = np.linspace(0, 2*np.pi, a)
    return [x0 + r*np.cos(t), y0 + r*np.sin(t), z0*np.ones_like(y0 + r*np.sin(t))]

class Reimann(HasTraits):
    radius=Range(1,10,5)
    numSlice=Range(1,100,5)
    scene = Instance(MlabSceneModel, ())

    def __init__(self,**traits):
        HasTraits.__init__(self,**traits)
        self.circle()
        self.inscribed_line()
        self.inscribed_chord()
        self.suprimum_line()
        self.suprimum_chord()
        self.find_area()

    def circle(self):
        x, y, z=circlepoints(self.radius,self.numSlice)
        #draw circle
        self.scene.mlab.plot3d(x, y, z,tube_radius=None)

    def inscribed_line(self):
        i=0
        g=0
        x,y,z=circlepoints(self.radius,self.numSlice)
        while i < self.numSlice:
            self.scene.mlab.plot3d([x0,x[g]],[y0,y[g]],[z0,z[g]],tube_radius=None)
            i=i+1
            g=g+100

    def inscribed_chord(self):
        i=0
        g=0
        b=0
        x,y,z=circlepoints(self.radius,self.numSlice)
        while i < self.numSlice:
            self.scene.mlab.plot3d([x[b],x[g]],[y[b],y[g]],[z[b],z[g]],tube_radius=None)
            i=i+1
            b=g
            g=g+100
        self.scene.mlab.plot3d([x[b],x[0]],[y[b],y[0]],[z[b],z[0]],tube_radius=None)
        #area of inscribed triangle
        base=np.sqrt((x[b]-x[0])*(x[b]-x[0])+(y[b]-y[0])*(y[b]-y[0])+(z[b]-z[0])*(z[b]-z[0]))
        self.area1=(base*self.radius)/2
        self.total_inscribed_area=self.numSlice*self.area1

    def suprimum_line(self):
        self.out_rad=self.radius/np.cos(np.pi/self.numSlice)
        i=0
        g=0
        x,y,z=circlepoints(self.out_rad,self.numSlice)
        while i < self.numSlice:
            self.scene.mlab.plot3d([x0,x[g]],[y0,y[g]],[z0,z[g]],tube_radius=None)
            i=i+1
            g=g+100

    def suprimum_chord(self):
        i=0
        g=0
        b=0
        x,y,z=circlepoints(self.out_rad,self.numSlice)
        while i < self.numSlice:
            self.scene.mlab.plot3d([x[b],x[g]],[y[b],y[g]],[z[b],z[g]],tube_radius=None)
            i=i+1
            b=g
            g=g+100
        self.scene.mlab.plot3d([x[b],x[0]],[y[b],y[0]],[z[b],z[0]],tube_radius=None)
        #area of suprimum triangle
        base=np.sqrt((x[b]-x[0])*(x[b]-x[0])+(y[b]-y[0])*(y[b]-y[0])+(z[b]-z[0])*(z[b]-z[0]))
        self.area2=(base*self.radius)/2
        self.total_suprimum_area=self.numSlice*self.area2

    def find_area(self):
        #area of original circle
        self.scene.mlab.text3d(-12,-13,-12,text="Area of infremum triangles  ="+str(self.total_inscribed_area))
        self.scene.mlab.text3d(-12,-15,-12,text="Area of suprimum triangles  ="+str(self.total_suprimum_area))
        average=(self.total_inscribed_area+self.total_suprimum_area)/2
        #print "Area of original circle"
        self.scene.mlab.text3d(-12,-17,-12,text="Both average ="+str(average))
        circle_area_formula=np.pi*self.radius*self.radius
        self.scene.mlab.text3d(-12,-19,-12,text="Area of circle using formula pi*r*r ="+str(circle_area_formula))
        #finding the area of one chord
        angle=(2*np.pi)/self.numSlice
        chord_area=((self.radius*self.radius)*(angle-np.sin(angle)))/2

        #find the total area of slice of suprimum triangle which is not inside the circle
        q=self.area2-self.area1
        t=chord_area-q
        small_area=self.numSlice*t
        self.scene.mlab.text3d(-12,-21,-12,text="Area of small part of outer triangle which is not inside the circle="+str(small_area))

    @on_trait_change('radius,numSlice')
    def update_plot(self):
        self.scene.mlab.clf()
        self.circle()
        self.inscribed_line()
        self.inscribed_chord()
        self.suprimum_line()
        self.suprimum_chord()
        self.find_area()

    # the layout of the dialog created
    view = View(Item('scene', editor=SceneEditor(scene_class=MayaviScene),
                    height=600, width=800, show_label=False),
                HGroup(
                        '_', 'radius', 'numSlice',
                    ),
                )

reimann=Reimann()
reimann.configure_traits()

this will produce various outputs depends upon changing radius and number of slice value

with text of outputs

This script have some camera object problem. It will clear soon. I want to say thanks to Enthought developers who helped me to do this and my guide Mr. Thiyagarajan.

3 Comments leave one →
  1. January 5, 2011 5:49 am

    you rocking man,,,, you digging in mayavi script and making gui using traits script…

    Keep it up…

    We expecting more from you man…

  2. raj kumar permalink
    February 8, 2011 2:27 pm

    keep research

  3. logeshprabu permalink
    February 8, 2011 2:38 pm

    Too hard to work like u!

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s

%d bloggers like this: