Take the 2-minute tour ×
Mathematica Stack Exchange is a question and answer site for users of Mathematica. It's 100% free, no registration required.

First, a bit of a long introduction to my problem:

I only have a few weeks of Mathematica experience. I am creating a mathematica application that calculates some material properties of steel based on metallurgical functions that include multiple chemical elements. I currently have the chemistry stored in a vector, and the list of elements stored in a second vector of the same length. Here's an example:

elements={"C","Mn","S","P","Si","Cu","Ni","Cr","V","Cb","Mo","Sn","Al","Ca","B","Ti","N","O","Ce","Pb","Bi","W","As","Zn","Sb"};
chemistry={0.032`,1.2`,0.0259`,0.006`,0.04`,0.15`,0.06`,0.04`,0.001`,0.051`,0.167`,0.007`,0.035`,0.0001`,0.`,0.001`,0.0087`,0.0176`,0.0002`,0.0001`,0.`,0.`,0.004`,0.036`,0.`};

I created a function called elem to look up the element name in elements, and find the corresponding entry in chemistry (or any other chemistry list of proper format):

elem[chem_,element_] := chem[[Flatten[Position[elements,element]]]][[1]];

For example, if I want to look up Nickel ("Ni") from the chemistry list, I use the following to get 0.06:

elem[chemistry, "Ni"]

The neat thing about this is when you write a formula using the elem function, the list elements can be used to print the formula in a more user-friendly way for verification. For example, the Titanium to Nitrogen ratio is one of the simpler formulas of interest:

TiByN[chem_] := elem[chem, "Ti"]/elem[chem,"N"]

When called with chemistry, it returns a number (0.114943), but when you call TiByN with elements, it returns the formula Ti/N. This way, I can print a list of formulas in a format that the metallurgists can read when they audit my application.

My question is one of performance. I need to do an optimization involving these formulas, and I think this method of indexing the chemistry is one of my big bottlenecks. Is there any way to precompile a function like elem for the list elements so that the position doesn't have to be looked up for every element in every formula on every iteration of the optimization?

Thanks

share|improve this question
1  
see "Functions That Remember Values They Have Found" in the manual –  george2079 Apr 3 '13 at 15:32
1  
This thread is very closely related--the principal (only?) difference is that other objects ("integer partitions") are the arguments rather than strings--mathematica.stackexchange.com/questions/22309/…. –  whuber Apr 3 '13 at 15:46
1  
As you mentioned that you are relatively new to Mathematica I wonder whether you have any indication that more common ways to handle this won't do? One straightforward way is (close to Istvans answer): define your formulas e.g. as TiByN="Ti"/"N" and then insert numeric values with something like TiByN /. parameters where parameters = Thread[elements->chemistry] (you could use Dispatch for optimization). With such an approach you can very easily control whether and when to insert values. What you do looks a little like working against the (strong parts of the) system to me... –  Albert Retey Apr 3 at 16:51
add comment

4 Answers

up vote 17 down vote accepted

Simple version using a variant of memoization

While part of the answer I was going to give was already posted by Istvan, I will still post mine since the self-precomputing part was not part of Istvan's answer. The following will use the variant of memoization to precompute the dispatch table:

ClearAll[elem];
elem[chem_, element_] :=
  With[{dispatchTable = Dispatch[Thread[ elements -> chem]]},
    elem[chem, el_] := el /. dispatchTable;
    elem[chem, element]];

What happens here is that the general new definition is created upon the first call with a given property, and then that more specific definition is used:

elem[chemistry, "C"]

(*  0.032  *)

Making a function generator

One flaw of the above answer is that the function elem depends on a global variable elements, in a rather non-transparent way. A better way would be to make a function generator:

ClearAll[makeElem];
makeElem[elemName_Symbol, elements_] :=
  elemName[chem_, element_] :=
    With[{dispatchTable = Dispatch[Thread[ elements -> chem]]},
      elemName[chem, el_] := el /. dispatchTable;
      elemName[chem, element]];

In this case, everything is local. One can now generate the elem function as

Clear[elem];
makeElem[elem, elements]

and then use it as before:

elem[chemistry, "C"]

(*  0.032  *)

Making this object-oriented

While the function generator is a better answer than the simple version, it may have a problem in that you have to come up with a name of the function elem - while this function is intimately tied with the particular set of elements. This calls for an OO approach, where it will be encapsulated in an instance (object), which will be created from a given list of elements. Here I will show how one can do this using an object-oriented extension I have proposed some time ago.

First, to try it, you have to execute

Import["https://gist.github.com/lshifr/4266126/raw/OO.m"]

which will load the OO package. Here is then a possible way to do what we want:

DeclareType[ChemicalElement][
   OO`Methods`new[elements_List] :=
     Module[{elem},
        $self@AddMethods[
              OO`Methods`addProperty[prop_] :=
                 Module[{dispatchTable = Dispatch[Thread[ elements -> prop]]},
                    $self@AddMethods[
                   OO`Methods`get[prop, el_] := el /. dispatchTable
                ]
             ]
        ];
        $self
     ]
];

It is perhaps a bit more complex than it needs be, since some of the mechanisms of the current OO` implementation are not perfect. What is done here is that I declared a type ChemicalElement with a new method. This method, when called, constructs an instance of ChemicalElement type and adds a new method addProperty (this is possible because the object model of OO` is flexible enough to be able to add methods on per-instance basis at run-time - which is e.g. not possible in Java but possible in Javascript).

Whenever you want to add a new property to your object, you call addProperty. So, here is our example:

ch = ChemicalElement[{}] @ new[elements]

(* « ChemicalElement[]» *)

Now we add the property chemistry:

ch @ addProperty[chemistry]

and now we can use it:

ch@get[chemistry,"C"]

(* 0.032 *)

If you then need to add another property (for example, and index of a given element), you can do it:

index = Range[Length[elements]];
ch@addProperty[index]

and now you can use it too:

ch@get[index,"C"]

(* 1 *)

One remaining flaw is that the properties are associated with specific values of specific variables (chemistry and index here). This can be easily changed by using some e.g. string names and changing the signature of addProperty accordingly.

Permanent installation of OO.m.

To permanently install OO` on your machine, you can use

Import["https://raw.github.com/lshifr/ProjectInstaller/master/BootstrapInstall.m"] 

which will install the ProjectInstaller project into your $USerBaseDirectory/Applications, and then call

Needs["ProjectInstaller`"] 
ProjectInstall[URL["https://gist.github.com/4266126/download"]]

alternatively, you can use the copyRemote utility by Rolf Mertig. Note however that ProjectInstaller installs entire project, which includes not just OO.m, but also example notebook and a file with project meta-data, so I recommend using that.

share|improve this answer
    
Thanks Leonid, this is exactly what I was looking for. –  ahhhh Apr 3 '13 at 16:14
    
@ahhhh Glad I could help, and thanks for the accept. See also my edit with an object-oriented approach. –  Leonid Shifrin Apr 3 '13 at 16:18
add comment

Use a dispatch table. It is an optimized element -> value lookup table that can be used to replace an element any time with its value. Now it does matching-and-finding every time, but if your list is not too big, this is pretty fast.

dispatch = Dispatch[Thread[elements -> chemistry]];

ratio[elemA_, elemB_, disp_] := (elemA/elemB) /. disp;
ratio[elemA_, elemB_] := ratio[elemA, elemB, {}];

ratio["Ti", "N", dispatch]
0.114943
ratio["Ti", "N"]
"Ti"/"N"
share|improve this answer
    
Thanks István, your answer was very educational as well, but I chose Leonid's because his was a drop-in replacement for my elem function. –  ahhhh Apr 3 '13 at 16:12
    
@ahhh No problem, and welcome to the site! –  István Zachar Apr 3 '13 at 16:38
add comment

You can use indexed variables with a shortcut function to define a long list of indices and values, like so:

indexedVariable[var_, indicesValues__] := With[{dict = List[indicesValues]}, Scan[(var[#[[1]]] = #[[2]]) &, dict]];
indexedVariable[chemistry, {"C", 0.032`}, {"Mn", 1.2`}, {"S", 0.0259`}];
chemistry["C"];

This will let the engine handle optimization of your lookups as it would any other indexed variable.

This method gleaned from http://www.physic.ut.ee/~kkannike/english/prog/mathematica/ under "TRUE DICTIONARIES IN MATHEMATICA 5", which includes many other helper functions for getting lists of keys and values and pairs which could be used to replicate your functionality of printing formulas numerically or symbolically.

share|improve this answer
add comment

With the new Association data structure introduced in the Wolfram Language/Mathematica 10 (you can try it now on the Raspberry Pi), this becomes extremely very simple to write and lookups are highly efficient as well.

property = Association@Thread[elements -> chemistry]
property @> Ni
(* 0.06 *)
share|improve this answer
    
Interesting but I'm wondering whether there are any advantages over simply chem=Thread@[elements-> chemistry]; Ni/.chem ? –  David Carraher Apr 3 at 16:23
    
@DavidCarraher With Dispatch, it comes very close to Association's timings. However, I like the cleanliness of using Association, especially with how the code reads. See also this answer and Leonid's comments under it for the advantages of Association. –  rm -rf Apr 3 at 16:34
    
Very useful link. Thanks. –  David Carraher Apr 3 at 16:50
add comment

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.