Astronomy 421 – Assignment 2

Name _____________________________________________ Date ______________

Astronomy 421 IDL Assignment 2

Writing procedures and functions

We are going to work through problems 3.9 and 3.12 of Carroll and Ostlie (2nd ed) through programming of the equations. As you follow these directions, pay close attention to how we program the problems and how they are solved. Please answer all questions posed within these directions.

First, some programming-style hints:

don't forget that IDL array indices start at 0
split programs into separate procedures and functions; test each separately making sure it works
use parameters instead of numbers to get a more easily adapted program
use SIZE(array) to get array dimensions for general procedures
choose clear variable names (IMPORTANT!)
add lots of comments (in ENGLISH please)
add detailed explanation at procedure/subroutine start
indent BEGIN ... END structures
JOURNAL, 'filename' records all your command-line entries, useful for subsequent conversion into programs

WHAT TO TURN IN (or, preferably, emailed to me): Copies of the programs you put together for steps I - V, if they differed from what is given in this assignment. Answer any questions posed below within the body of the email or within the context of the programs themselves. There should be 3 copies of chapters 2, 3, and 4 of Liam Gumley's Practical IDL Programming in the computer lab room to help you will all of this, and to lead you to even greater accomplishments in IDL.

Carroll & Ostlie, Problem 3.9

Step I: Using a text editor ( not MS Word or similar formatting editor, see Linux tutorial), edit a file called problem3_9.pro and copy and paste in the following C & O problem3_9.pro and save the file to the directory of choice. You may also wish to copy my file that contains a short list of astronomical physical constants: physical_constants.inc. If you use this file, then you will need to remember to use those specific variable names.

NOTE: The command help used in this IDL program is a way to get all of the variables printed out, and is good for debugging, but not for printing of specific variable values.

Step II: Now, start IDL and compile the program within IDL by typing

IDL> .compile problem3_9 ; IDL assumes a .pro suffix

IDL> problem3_9 ; tell IDL to run it! This is how we call procedures and functions

You should always work through part of a problem the long way (using a calculator) to recognize that is always an essential task to check the calculations that are coming out of your computer! If your program is solving the same type of problems repeatedly, then you need only to spot check. But, you must never simply accept what your computer produces. How will you know if the answer is right? Did your long hand answers agree with your computer answers?

Step III: Continue with good programming techniques and write a NEW program to find out the peak wavelength, in nanometers, for a given temperature using Wien's Law where you are prompted to enter the value of the temperature from the terminal. (In other words, don't "hard-wire" the value in the program. Be sure to save the commands. You can do this as a command-line process using the JOURNAL command.

Step IV: Continue with good programming techniques and write a program that returns the value of the Planck function for just one temperature (of your choosing) at a wavelength of 500 nanometers. You can again do this via command lines and JOURNAL to test values and results as you go. Test that your answer is correct by using your calculator or checking your answer with a classmate who programmed independently from you.

Step V: Now, we are going to incorporate the Wien program as a FUNCTION and Planck program as a PROCEDURE into the main program. Name your Wien program wien.pro and the Planck program findplanck.pro [WHY can't we call the procedure 'planck.pro'?). Your modified functions will need to use the values that are set in the main procedure, so get rid of any prompting for entering values via the terminal.

Edit these programs so that the first and last lines read as follows (where you insert your commands where I have typed ;statements ):

FUNCTION WIEN, TEFF_STAR
; statements
END
;
;

PRO FINDPLANCK, TEFF_STAR, PLANCK_VALUE
; statements
END

Edit your main program, problem3_9.pro , by putting in calls to the function and procedure. You do that by adding the following lines before the print statement (don't forget to print out these new variable values):

wave_max = wien(teff_star) ; parentheses indicate parameters passed to function

findplanck, teff_star, planck_value

You must now recompile and run the main procedure. Do this in IDL. When you .COMPILE you must include all of the file names that are called. You can also compile them separately. [NOTE: In past quarters, the version of IDL we were running insisted that any procedures or functions called by the main procedure had to be compiled first.]

You will most likely get errors! Try a few minutes to solve them. Read through Ch. 3 Writing IDL Programs handout. Don't break bricks with your head if you can't figure it out have a classmate help you or your instructor. Email works if it is outside of class time. Do you have undefined variables? For now, you must redefine variables within each procedure or function, since they don't automatically know the variables. There is a better way; see Ch. 3 of Gumley (p. 110 - 128 ).

Carroll & Ostlie (2nd ed), Problem 3.12: Solution

What to turn in: 1) a summary of your comparison of the programs, 2) answers to the questions at the end (especially the question that asks, "What have you found?").

An extremum in the Planck function (Eq. 3.22), occurs when the derivative is equal to zero: .

We have to differentiate by parts:

which simplifies to

Defining , this becomes , which cannot be solved analytically. We must solve numerically, finding an x such that .

This solution requires a root-finding program . I chose to use two of the most simple ones that are part of the Press et al Numerical Recipes in Fortran 77. It is quite easy to convert Fortran code into IDL.

I found that the old, slow work horse rtbis.f (now, rtbisidl.pro ) that finds a root via the bisecting method (you must know the range in which the root lies) gave the correct answer with the least amount of fussing. Come to find out, the equation for which I needed to find the root actually has two roots, one of which is 0 . The second root-finding program (the one I usually prefer as it is faster for complicated equations) insisted on converging to the x = 0 root it got away from the correct root and just kept on going. Once I bracketed the known, non-zero root more tightly, the program worked great. You should ALWAYS be aware that a program may go off into territory that is uncharted, and return answers that you do not want. The second program is from rtsec.f (now rtsecidl.pro ), and uses a secant method to converge on a solution.

Copy the following files to your IDL directory:

/net/projects/Astro_421/IDL/rtbisidl.pro

/net/projects/Astro_421/IDL/rtsecidl.pro

You review the programs and try to figure out the logic that is used, as well as the commands. You will need all those IF THEN BEGIN ENDIF ELSE ENDELSE ENDFOR's later! You should copy them into your program directory, compile each one in turn and run it.

How do the programs compare?

Which one converged on a solution with the least number of iterations?

Your solution gives you the value for . Solve for the wavelength in terms of the inverse temperature, substituting the values for the constants. Be sure to put down the correct units. What have you found? Did you do all of this within IDL, the easiest way?