Why is Mathematica Stopping My Long Calculation?
January 30 2016 in Mathematica, Programming | Tags: | Author: Christopher Rackauckas
You set it all up, you know it’s going to take a few days, and you run it. You check it a few hours later, it’s all good! The next morning… there are no cells computing still?
Mathematica’s documentation isn’t that helpful for what just happened. If you were using a function like FullSimplify, there are time limits on it. However, it will give you an error/warning if it hits the time limit, so if you ended up with a blank screen with no calculations, that’s not it.
This happened to me. I found out that the hidden culprit (for me) was Mathematica’s history tracking. The fix is simple, add the following code to the top of your file:
$HistoryLength = 10;
What’s happening is that Mathematica saves its entire command history by … READ MORE
Julia iFEM3: Solving the Poisson Equation via MATLAB Interfacing
January 24 2016 in FEM, Julia, MATLAB, Programming | Tags: FEM, julia, Poisson Equation | Author: Christopher Rackauckas
This is the third part in the series for building a finite element method solver in Julia. Last time we used our mesh generation tools to assemble the stiffness matrix. The details for what will be outlined here can be found in this document. I do not want to dwell too much on the actual code details since they are quite nicely spelled out there, so instead I will focus on the porting of the code. The full code is at the bottom for reference.
The Ups, Downs, and Remedies to Math in Julia
At this point I have been coding in Julia for over a week and have been loving it. I come into each new function knowing that if I just change the array dereferencing from () to [] and wrap vec() calls around vectors being used as … READ MORE
Julia iFEM 2: Optimizing Stiffness Matrix Assembly
January 23 2016 in FEM, Julia, MATLAB, Programming | Tags: julia, optimization, sparse, vectorization | Author: Christopher Rackauckas
This is the second post looking at building a finite element method solver in Julia. The first post was about mesh generation and language bindings. In this post we are going to focus on performance. We start with the command from the previous post:
node,elem = squaremesh([0 1 0 1],.01)
which generates an array elem where each row holds the reference indices to the 3 points which form a triangle (element). The actual locations of these points are in the array node, and so node(1) gives the points in the (x,y)-plane for the $$i$$th point. What the call is saying is that these are generated for the unit square with mesh-size .01, meaning we have 10201 triangles.
The approach to building the stiffness matrix for the … READ MORE
Optimizing .*: Details of Vectorization and Metaprogramming
January 21 2016 in Julia, MATLAB | Tags: BLAS, de-vectorization, high performance computing, Linpack, MKL, VML | Author: Christopher Rackauckas
For many of us mathematicians, we were taught to use MATLAB, and we were taught to vectorize everything. I mean obviously if we have matrices $$A$$, $$B$$, and $$C$$ and want to multiply element-wise (say to solve a reaction-equation at each point in space), then the optimal code is
A.*B.*C
No questions to ask, right? Actually, this code isn’t as optimized as you’d think. Lets dig deeper.
BLAS, Linpack, and MKL
The reason you are always told by “the lords of numerical math” to vectorize your code is because very smart programmers worked really hard on making basic things work well. Most of the “standardized” vectorized computations are calling subroutines from packages known as BLAS and LINPACK. To see what version your MATLAB is using, you can call
Quick Optimizations in Julia for Performance: A Practical Example
January 19 2016 in Julia, MATLAB | Tags: AVX512, julia, performance, SIMD, threading | Author: Christopher Rackauckas
Let’s take a program which plots the standard logistic map:
r = 2.9:.00005:4; numAttract = 100; steady = ones(length(r),1)*.25; for i=1:300 ## Get to steady state steady = r.*steady.*(1-steady); end x = zeros(length(steady),numAttract); x[:,1] = steady; for i=2:numAttract ## Now grab some values at the attractor x[:,i] = r.*x[:,i-1].*(1-x[:,i-1]); end using PyPlot; fig = figure(figsize=(20,10)); plot(collect(r),x,"b.",markersize=.06) savefig("plot.png",dpi=300);
This plots the logistic map. If you take the same code and change the array … READ MORE
Julia iFEM1: Porting Mesh Generation
January 19 2016 in FEM, Julia, MATLAB, Programming | Tags: iFEM, julia, MATLAB, matplotlib, triangulation | Author: Christopher Rackauckas
My first project on the quest for a Julia finite element method is a simple homework problem. Just some background, this is for UC Irvine’s graduate Computational PDEs 226B course where in the first quarter we did all sorts of finite difference methods and now is our first foray into finite element methods. The purpose of the project is to grasp the data structure enough to use simple tools (i.e. mesh creation and plotting) to create a finite element solver for Poisson’s equation in 2D and check the performance differences. For those who haven’t programmed much this is a great learning exercise, but being a pretty standard exercise the MATLAB code took no time to writeup and so I started thinking: how does Julia compare to MATLAB for solving simple PDEs with the finite element method?
Testing this would take … READ MORE
WordPress Julia Syntax Highlighting
January 19 2016 in Julia | Tags: GeSHi, julia, wordpress | Author: Christopher Rackauckas
using Julia print("Test for Julia") A = [0 1 ; 0 1];
Syntax highlighting in Julia exists? Yes, you can make it work even though no plugins currently support it. The way I was able to get it to work was to use the WP-GeSHi-Highlight WordPress plugin which uses GeSHi to perform the syntax highlighting. Sure enough, an unofficial Julia syntax file for GeSHi is hosted on GitHub. I manually plopped this file into the GeSHi folder in my WordPress install, used the command lang=”julia” and you can see at the top of this page that it worked!