Sunday 18 August 2013

GSoC Week9 : Beginning of electrostatics

Week 8 of GSoC is now over, and just 5 more weeks remain till the hard pencils-down date.
However, my work on the project will most likely exceed that time, since I will be having my first tests around September 17. I am sure I will be able to get all the code in by that time, but to make the entire work shippable, things like Sphinx documentation, doctests etc will also have to be finished, which may take a little more time beyond the actual deadline.

Anyways, coming to this week's work, I wrote some code for electrostatics. I added some functions and the ParticleCharge class (inheriting from Particle of mechanics module).
However, the E-M module's handling of charged particles required me to modify my earlier implementation of Particle- most importantly, adding a 'set_motion' method. Like I mentioned earlier, Particle earlier had a frame attached to it. Gilbert and I discussed at length on this, and we agreed that making the user initialize a new frame everytime, just to attach it to a Particle, was quite cumbersome.
Hence, I modified the class in the mechanics core- It still has a frame as an attribute, but its now 'hidden' from the user (nothing in Python in hidden from the user, but what I mean here is that the user can do everything with a Particle without caring about the frame). The set_motion method I mentioned earlier just initializes a new frame according to the motion parameters set by the user (similar to MovingRefFrame).

The main reason I want a frame attached to a Particle is the ease of inner working. Since a Particle will not have an 'angular acceleration', those parts of the MovingRefFrame code will never really get executed in Particle. However, Gilbert has a point in pointing out that its unnecessarily adding a complex attribute to the class. We havent finalized anything yet, but re-writing the motion-setting methods for Particle would be tedious, and more importantly, will lead to quite a lot of code duplication. Lets see what we end up doing.

Anyways, I also finished writing tests for the code I wrote this week, and I am waiting for a review from Gilbert. Jason and he also encouraged me to add example-tests to ensure everything works smoothly.
Here are the example tests I wrote (they are really basic, but its better that way, since strenuous tests only make things worse)-

 def test_example_1():  
   """  
   This is a sample example meant to test the basic functionality of  
   electrostatics.py  
   
   Consider a frame R, and a particle charge P(mass m, charge q) situated  
   at its origin initially.  
   Now, at time t=0, a field R.x+R.y+R.z is 'switched on' and kept that way  
   for time 't0'.  
   The user needs the x-coordinate of the particle charge as a function of  
   time, after the field is switched off at time=t0.  
   """  
     
   #Basic initialization  
   m, q = symbols('m q')  
   P = ParticleCharge('P', m, q)  
   P.set_motion(R, pos_vector = 0)  
   field = R.x + R.y + R.z  
   time = dynamicsymbols._t  
   t0 = Symbol('t0')  
   #The acceleration is equal to the electrostatic force experience by the  
   #particle charge, divided by its mass  
   acceleration = P.electrostatic_force(field)/P.mass  
   #Use get_motion_acc from the mechanics core to find velocity and position  
   #parameters using acceleration function and boundary conditions  
   translation = get_motion_acc(acceleration, \  
                  0, P.pos_vector_wrt(R), frame = R)  
   #Calculate the motion parameters of the particle charge at time t0  
   acceleration = translation[0]  
   velocity_at_t0 = translation[1].subs({time:t0})  
   position_at_t0 = translation[2].subs({time:t0})  
   #Set motion of P accordingly  
   P.set_motion(R, trans_acc = 0, trans_vel_b = velocity_at_t0,  
          pos_vector_b = position_at_t0)  
   #assert that the dot product of P's pos_vector wrt R, with R.x  
   #matches the actual result  
   assert P.pos_vector_wrt(R).dot(R.x) == \  
       q*t*t0/m + q*t0**2/(2*m)  
   
   
 def test_example_2():  
   """  
   This is a sample example meant to test the basic functionality of  
   electrostatics.py  
   
   Consider a set of 3 particle charges of equal mass and charge placed  
   symmetrically on a circle of unit radius around the origin, with one  
   point lying at position vector of R.y  
   We calculate the energy required to assemble this system of charges from  
   infinity(with no kinetic energy in any particle) and find out values of  
   the electrostatic fields and potentials at the origin  
   """  
   
   #Basic initialization  
   from sympy import sin, cos, pi  
   m, q = symbols('m q')  
   P = ParticleCharge('P', m, q)  
   Q = ParticleCharge('Q', m, q)  
   S = ParticleCharge('S', m, q)  
   #Place particle charges as per required configuration  
   P.set_motion(R, pos_vector = R.y)  
   Q.set_motion(R, pos_vector = -cos(pi/6) * R.x + -sin(pi/6) * R.y)  
   S.set_motion(R, pos_vector = cos(pi/6) * R.x + -sin(pi/6) * R.y)  
   #Check outputs of required values  
   assert charge_assembly_energy(P, Q, S) == 3*3**(-0.5)*k*q**2  
   assert S.electrostatic_field(R, 0) + \  
       Q.electrostatic_field(R, 0) + \  
       P.electrostatic_field(R, 0) == 0  
   assert P.electrostatic_potential(R, 0) == \  
       Q.electrostatic_potential(R, 0) == \  
       S.electrostatic_potential(R, 0)  

If you think the API is a little complex, I think so too, but those modifications are for later. It would be good, according to me, to have something like an 'environment' infrastructure. Basically, the user will keep adding components (fields, particles etc etc) to it, and the constituents will keep getting modified accordingly. However, it will be good to have the module at a certain level of working before thinking on those lines.

That's all there is to report for this week. The coming week, I will be reading up on magnetostatics from Griffith's and start coding for the same. However, I am afraid things are only going to get more and more complex hereon. Lets hope for the best.

Have a great week :-)

Sunday 11 August 2013

GSoC: Week 8

Week 8's over.
This week was more of 'studying' and 'preparing' rather than coding. Since the entire infrastructure behind my initial proposal has now shifted considerably, I am drawing out the plans anew and getting the basic things done before I start working on sympy.physics.em (tentative name).

First off, I finally got the logic docs up at the SymPy development docs. You can see them here.
Second, I wrote the basic field functions for electrostatic theory (including the tests).
Apart from that, I was mostly busy with college stuff (first week of classes are usually a little difficult for adjusting) and going through Griffith's book. So hopefully, I will get a lot done by the end of the coming week (code + further plans).
Prasoon's code is shaping up well, and I may help him with the testing and debugging part to speed up the development of sympy.vector. I feel quite optimistic about the projects (Prasoon's and mine) now.

Well, thats all (short post, I know) for now.
Will come back with more next time.
Have a great week :-)

Sunday 4 August 2013

GSoC Midterm Report (Week 7)

Week 7 of GSoC is over, which means I am now in the 'second half' of the coding period. The midterms were conducted by Google from 29th to 2nd, and I am happy to report I passed them :-).

As I had predicted in my last blog post, this week wasn't very exciting, atleast as far as the project is concerned. I am now at Goa, and classes have started. However, I did get a few things done - The Dyadic framework is now almost complete, with all the classes functioning as expected (from point of view of dyadic operations). I also added all the special operations like cross, dot, express and time-differentiation. For now, just like last week, they use SymPy Symbols for temporary working.
Moreover, I finished writing all the tests for my first PR. Gilbert reviewed them, and I made a few changes (mostly API-wise).
I also tested my MovingRefFrame implementation of vector time-derivatives on my hacked mechanics module, and I am happy that it worked fine. Thus, I now have the proof-of-concept for all the work I have done so far.
I have also submitted a PR to get the logic module docs into the online development documentation, but I am quite sure I am going wrong somewhere. I have pinged Aaron, and I hope it gets solved by the end of this week.

To sum up my work till now, this is what I got done-

1) MovingRefFrame class - The basic class to represent moving reference frames in 3-D space. Inherits from the CoordSysRect class of the new vector module.

2) Particle class - Class to represent particles (mass but no volume) in space. Each particle is associated with a frame of its own, and all it's parameters are calculated wrt it.

3) Miscellaneous mechanics functions like the ones to obtain position vector given a time-dependent vector and boundary conditions.

4) Dyadic framework - Classes Dyadic, BaseDyadic, DyadicAdd, DyadicMul. All these classes are based on relevant SymPy core classes and build dyadic-based functionality on top of them.
Currently, I am reading Griffith's book (very helpful suggestion from Aditya - my collegemate) and taking notes for the EM module. Since the development of the new mechanics core is now nearing completion, I can focus time on working on the new addition to sympy.physics.

Hence, the coming week will be spent on polishing off the mechanics work, making sure all required parts are in place, maybe helping Prasoon with some functions of vector fields (so that he can concentrate on vector-space integration).
Plus, I have to design the vector-dyadic interfacing. Basically, how dyadics will be instantiated from vector 'outer' functions. Obviously, having a dyadic framework is useless without fitting it to the vector module.
And last but not the least, I plan to code the electrostatic functions for the EM sub-module by this weekend. Hope I succeed.

Have a great week :-)