# Wanted to write a mathematica program

1. Jun 21, 2005

### steven187

hello all

well I wanted to write a mathematica program so i can be able to experiment but for some reason i cant get it to work, does anybody have any ideas on where im going wrong, below is something that i have cut out of my mathematica program, all you have to do is cut and paste into mathematica

thank you

\!$$\(LowerRiemann[a0_, b0_, n0_] := \[IndentingNewLine]Module[{a = a0, b = b0, n = n0, X}, \[IndentingNewLine]\[IndentingNewLine]X\_i_ = a + i\ \(b - a$$\/n\ ; \ \[IndentingNewLine]sum =
0; \[IndentingNewLine]For[\[IndentingNewLine]i = 1,
i <= n, $$i++$$, \[IndentingNewLine]Xv =
Select[Xvals,
X\_$$i - 1$$ <= #1 <= X\_i\ &]; \[IndentingNewLine]Ymin = \
Min[f[Xv]]; \ \[IndentingNewLine]sum =
sum + Ymin\ $$b - a$$\/n;\ \ \[IndentingNewLine]]; \ \
\[IndentingNewLine]Return[\ \ sum\ ];\ ];\)\ \[IndentingNewLine]
$$f[x_] = x\^2;$$\[IndentingNewLine]
$$a = 0;$$\[IndentingNewLine]
b = 1; \[IndentingNewLine]
$$n = 25;$$\[IndentingNewLine]
LowerRiemann[a, b, n]\)

2. Jun 22, 2005

### steven187

hello all

LowerRiemann[a0_, b0_, n0_] :=

Module[{a = a0, b = b0, n = n0, X},
X_k= a + k (b - a)/n ;
sum = 0;
For[
i = 1, i<=n, i++,
Xv = Select[Xvals, X_{i-1}<=#1<=X_i &];
Ymin = Min[f[Xv]];
sum = sum + Ymin( b - a)/n;
];
Return[ sum ];
];

f[x_] := x^2;

a = 0;
b = 1;
n = 25;

LowerRiemann[a, b, n]

3. Jun 22, 2005

### saltydog

Jesus Steven. I though I was good with Mathematica! I do not see where you're defining the array Xvals. And would you do me a big favor? Could you please explain the # and & and in general what's going on in this:

Xv = Select[Xvals, X_{i-1}<=#1<=X_i &];

4. Jun 22, 2005

### saltydog

This is how I would calcualate a lower Reimann sum:

LowerReimann[xstart_,xend_,f_,n_]:=Module[{i,xinc,xdif,xcurrent,sum},

xinc=xstart;
xdif=(xend-xstart)/n;
xcurrent=xstart-xdif;
sum=0;

For[i=1,i<=n,i++,
xcurrent+=xdif;
yval=FindMinimum[f[x],{x,xcurrent,xcurrent,xcurrent+xdif}];
sum+=yval*xdif;
];
Return[sum];
];

Last edited: Jun 22, 2005
5. Jun 22, 2005

### steven187

hello all

well saltydog about this Xv = Select[Xvals, X_{i-1}<=#1<=X_i & ]; this is a logical arguement, the "Select" picks out elements out of the list of values "Xvals" which satisfy the criterion X_{i-1}<=#1<=X_i & where #1 stands for the first arguementand the & is the intersection of the arguement in which in this case are intervals and so this "Select[Xvals, X_{i-1}<=#1<=X_i & ];" is used by the Xv as a temporary variable to carry out its iterations, from what i remember I dont think Xvals needs to be declared it is not a local variable which is significant to the module but im not sure I learnt this stuff years ago, now in terms of your code it didnt work, im going to look into it, and i will update you if there are any changes, did you get your code to work?

steven

6. Jun 22, 2005

### saltydog

Alright Steven. You caught me. I didn't check it. Silly me. I should know better. I'll perfect it then post my version of the code. Really, once we get it perfected, I'll like to do the Upper limit and then run up the iterations real high and verfy numerically that they converge for Reimann-Integrable functions.

7. Jun 22, 2005

### steven187

well saltydog

im thinking that we shall join the upper and lower code and build a convergence tester such that it tests the error between the two as it should get smaller and smaller as we increase the number of partitions.any way lets get this thing workin

steven

8. Jun 22, 2005

### saltydog

Yea, right Steven. And I though I was a good programmer. But that's Ok. I am, and will figure it out. Castilla, are you reading this? I hope you're happy.

9. Jun 22, 2005

### steven187

well saltydog, you are a good programmer no doubt about that, your style of programming is nice and understandable, I had mastered C programming,but once i discovered how quik you can program things in mathematica it became my ideal program for numerical problems, the only difference is that mathematica has built in functions where C can tell you where you went wrong well my one does ;) by the way have i created an error anywhere with what i have mentioned?

10. Jun 22, 2005

### saltydog

Hello Steven. I cannot understand the syntax you're using with Select and Mathematica seems to have difficulty too with it. I'm quite stubborn and disappointed Mathematica does not have a simple means of determining the min or max of a function on an arbitrary interval. If you know, please let me know. I am however, a very practical person. In this regards, I wrote my own (inefficient) algorithm to determine such called GetFunctionValue. Supplying either "HIGH" or "LOW" gets the value in the indicated interval. The rest of the code used to check x^3 is below as well as a plot of the difference between the high and low values which of course will approach 0 as n is increased.

Code (Text):
(* "HIGH" flag returns high value, "LOW" returns low *)

GetFunctionValue[flag_, funct_, start_, end_] :=
Module[{xval, dif, factor, savevalue},
If[end < start,
Return[0];
];
savevalue = 0;
factor = 0.01;
dif = Abs[(end - start)] factor;
savevalue = funct[start];
For[xval = start, xval <= end, xval += dif,
If[flag == "HIGH",
If[funct[xval] > savevalue,
savevalue = funct[xval];
];
,
If[funct[xval] < savevalue,
savevalue = funct[xval];
];
];
];
Return[savevalue];
];

(* calculate upper or lower Reimann sums *)

Reimann[flag_, f_, xstart_, xend_, iter_] := Module[{i, deltax, sum, xval},
deltax = Abs[(xend - xstart)]/iter;
sum = 0;
xval = xstart;
For[i = 1, i <= iter, i++,
sum += (GetFunctionValue[flag, f, xval, xval + deltax]) deltax ;
xval += deltax;
];
Return[sum];
];

(* calculate difference between upper and lower sums *)

RDif[f_, xstart_, xend_, partnum_] :=
Reimann["HIGH", f, xstart, xend, partnum] -
Reimann["LOW", f, xstart, xend, partnum];

(* test case*)

ftest[x_] := x^3;
rvalues = Table[{n, RDif[ftest, 0, 5, n]}, {n, 1, 100}]

#### Attached Files:

• ###### reimanndif.JPG
File size:
4.5 KB
Views:
54
11. Jun 22, 2005

### steven187

hello saltydog

I kind of slept again while im at it :zzz: , luckily i did get sum where , and its because of your
doubts in defining Xvals, the thing is I forgot to include the partition and thats what Xvals is, now I will explain to you in the most simple terms about what select does in terms of this problem, ok select[list of elements,criterion] Xvals is a list of elements and it is also the list in order of all the endpoints of the subintervals, with this
X_{i-1}<=#1<=X_i & what this does is finds the elements which satisfy
X_{i-1}<=Xvals<=X_i, and there are only 2 elements X_{i-1} and X_i, and so for the first arguement "#1" when i equal to 1 we get x_0<=Xvals<=x_1, what the # does is provide a slot such that any value that satisfy this arguement will be assigned to Xv and so as a result we are not assigning one value for Xvals we are assigning a whole interval and thats why we use a slot # and the "&" character is there so that it groups all the points which satisfy the arguement and in this case it is a whole interval. if you have any further problems just update me, now in terms of your code il check it out, but when you said

disappointed Mathematica does not have a simple means of determining the min or max of a function on an arbitrary interval

thats what I have done with select it assignes Xv with an arbitrary interval,

now in terms of my code I added a line to define Xvals "thanxs Saltydog" in terms of a list of all the elements in the partition but I keep getting value 0 check it out

steven

LowerRiemann[a0_, b0_, n0_] :=

Module[{a = a0, b = b0, n = n0, X},
X_k= a + k (b - a)/n ;
sum = 0;
Xvals = Table[N[a + j(b - a)/n], {j, 0, n}];
For[
i = 1, i<=n, i++,
Xv = Select[Xvals, X_{i-1}<=#1<=X_i &];
Ymin = Min[f[Xv]];
sum = sum + Ymin( b - a)/n;
];
Return[ sum ];
];

f[x_] := x^2;

a = 0;
b = 1;
n = 25;

LowerRiemann[a, b, n]

Last edited: Jun 22, 2005
12. Jun 22, 2005

### saltydog

Hello Steven. I'm very much a purist when it comes to writing tite code and thus I prefer your "clean" method to my messy one. I'll spend time with it as I also want to improve my skills using Mathematica.

Also, if you want, you can just cut-and-paste Mathematica code into a post as long as it doesn't have "traditional formatted" codes. If it does, like exponents, fractions, just rewrite them in line-mode or else Mathematica will insert formatting code everywhere. Also, when you paste it to a post, wrap it with CODE tags (the # above), this will format it (indent) like it was in Mathematica. That's why my code above is in blue.

13. Jun 22, 2005

### saltydog

Hello Steven. I changed your code. I converted X to a function (maybe that's the same as you were doing), added the function as a parameter to the function call, made a few more changes. Seems to be working. Check it out. Let me know if not. The Select syntax and function is still very hard for me to fathom. I'll keep working on it because I think I'm a better programmer than Mathematician.

Also, slight discrepancy with the Table size: actually n+1 elements but only using n elements. Suppose we could clean up that too.

Code (Text):
LowerRiemann[func_, a_, b_, n_] := Module[{sum, Xvals, Xv, Ymin, i, k},
X[k_] := a + k (b - a)/n;
sum = 0;
Xvals = Table[N[a + j(b - a)/n], {j, 0, n}];
For[i = 1, i <= n, i++,
Xv = Select[Xvals, X[i - 1] <= #1 <= X[i] &];
Ymin = Min[f[Xv]];
sum = sum + Ymin (b - a)/n;
];
Return[sum];
];

a = 0;
b = 1;
n = 50;
func[x_] := x^2;
LowerRiemann[func, a, b, n]

14. Jun 23, 2005

### steven187

hello saltydog

well i have chucked your code into mathematica but i aint geting any numbers the only output i get is the sum of 50 functions, it looks like to me that your program occurs when it is selecting the 2 elements from the partition that satisfy the interval (which are the end points of each sub interval), and so it is suppose to be finding the minimum function values of an interval, but instead its trying to find the minimum of the function values of those 2 elements, but it wont calculate the function vales of those two points to find the minimum for some reason i think its got to do with how you have declared the function f(x), anyway im gonna try and fix it up and get it over and done with today hopefully, by the way im pretty sure that there is a difference between indexing and defining things as a function, one big difference is the domain that is involved, any way keep me updated

steven

15. Jun 23, 2005

### saltydog

Steven . . . I made a mistake with this line:

Ymin = Min[f[Xv]];

It should say:

Ymin=Min[func[Xv]];

Tell you what. I bet a dollar it works now.

16. Jun 23, 2005

### steven187

well good work saltydog i knew you would figure it out
well here is the upper :tongue2:
Code (Text):
UpperRiemann[func_, a_, b_, n_] := Module[{sum, Xvals, Xv, Ymax, i, k},
X[k_] := a + k (b - a)/n;
sum = 0;
Xvals = Table[N[a + j(b - a)/n], {j, 0, n}];
For[i = 1, i <= n, i++,
Xv = Select[Xvals, X[i - 1] <= #1 <= X[i] &];
Ymax = Max[func[Xv]];
sum = sum + Ymax (b - a)/n;
];
Return[sum];
];

a = 0;
b = 1;
n = 50;
func[x_] := x^2;
UpperRiemann[func, a, b, n]
Error= Abs[UpperRiemann[func, a, b, n] -LowerRiemann[func, a, b, n]]
which converges yay we did it

17. Jun 23, 2005

### steven187

well saltydog

I am now trying to manipulate this code such that it is capable of working with piecewise continuous or piecewise discontinuous functions such as
f(x)= 1 if x>=0
=-1 if x<0 by the way this is riemann integrable on the interval [-1,1] as from my understanding in another thread I had started, well Im gonna het to it il update you if I get anywhere

steven

18. Jun 23, 2005

### steven187

hello all

for some reason this code does not even work with piecewise continuous functions either, I have been working with this function

func[x_] := Which[x < 0, -1, x >= 0, x^2 - 1]

but it fails to give any answers only the sum of functions, have i defined this piecewise function incorrectly? any advice would be helpful

thank you