Advanced Search

Carl Love

Carl Love

28095 Reputation

25 Badges

13 years, 100 days
Himself
Wayland, Massachusetts, United States
Contact Carl Love
My name was formerly Carl Devore.

MaplePrimes Activity


These are answers submitted by Carl Love

More solutions...

Carl Love 28095
January 29 2014
0 6

Motivated by Axel's analysis, I looked more closely for solutions for the situations that RootFinding:-NextZero rejected. I used a "microscope", searching very close to the origin. By that I mean that I used semilogplot and applied fsolve to numer(f(10^d)) rather than to f(d). I plotted in all cases, and the sequence of plots really shows what's happening.

The first point remains a mystery because the limit(f(d), d=0, right) = -infinity (as Axel noted) where one would expect 0. The limit for the other points is in fact 0, and the plots seem like a continuous progression.

 

restart:

L := convert(
     [[1.0, 0.75e-1], [2.0, .1], [3.0, .12], [4.0, .14], [5.0, .16],
     [6.0, .18], [7.0, .2], [8.0, .22], [9.0, .24], [10.0, .26]],
     rational
):

NN := -N+(N0B-NB)*(erf((1/2)*x/(sqrt(t)*sqrt(d)))
     -sqrt(erf((1/2)*alpha/d)))/sqrt(erfc((1/2)*alpha/d))+NB:

 

alpha:= 2*10^(-12):  NB:= 75/1000:  N0B:= 2/10:  t:= 360000:

 

NNlog:= eval(NN, d= 10^d):

Digits:= 50:

for k to nops(L) do
     f:= numer(simplify(eval(NN, [x,N]=~ L[k])));
     f_log:= numer(simplify(eval(NNlog, [x,N]=~ L[k])));
     printf("\nTrying x= %f, N= %f\n", evalf[5](L[k])[]);
     dd:= fsolve(f_log, d= -6..-3); #10^(-6) .. 10^(-3)
     if dd::numeric then
          printf("Solution found: d= %a\n", evalf[15](10^dd));
     else
          printf("No solution found.\n");
     end if;
     print(plots:-semilogplot(f, d= 1e-14..1e-2, labels= [d,'NN']))
end do:


Trying x= 1.000000, N= 0.075000

No solution found.


Trying x= 2.000000, N= 0.100000
Solution found: d= .864547473017420e-4


Trying x= 3.000000, N= 0.120000
Solution found: d= .570969017295855e-4


Trying x= 4.000000, N= 0.140000
Solution found: d= .445134288859349e-4


Trying x= 5.000000, N= 0.160000
Solution found: d= .350843030752839e-4


Trying x= 6.000000, N= 0.180000
Solution found: d= .253007839385319e-4


Trying x= 7.000000, N= 0.200000

No solution found.


Trying x= 8.000000, N= 0.220000

No solution found.


Trying x= 9.000000, N= 0.240000

No solution found.


Trying x= 10.000000, N= 0.260000

No solution found.

 

Download erf_solve_2.mw

Law of cosines...

Carl Love 28095
January 29 2014
2 0


restart:

OA,OB,OC:= 3$3:  AB:= 3*sqrt(2):  BC:= 2*sqrt(2):

Law of cosines:

LoC:= c^2 = a^2 + b^2 - 2*a*b*cos(C):

eval(LoC, [c= AB, a= OA, b= OB, C= AOB]):

AOB:= solve(%, AOB);

(1/2)*Pi

eval(LoC, [c= BC, a= OB, b= OC, C= BOC]):

BOC:= solve(%, BOC);

arccos(5/9)

eval(LoC, [c= AC, a= OC, b= OA, C= 2*Pi-AOB-BOC]):

solve(%, AC);

14^(1/2)+2, -14^(1/2)-2

Take the positive value. So the answer is

%[1];

14^(1/2)+2

 


Download triangle.mw

limits of integration...

Carl Love 28095
January 26 2014
0 0

In order to get a numeric result, you need two numeric limits of integration. Use Int instead of Intat, and then apply evalf.

RootFinding:-NextZero...

Carl Love 28095
January 26 2014
2 1

RootFinding:-NextZero is not as strict as fsolve is with the accuracy. It obtains a solution for several of your (x, N) pairs. For the pairs for which no solution is found, a plot indicates that there is truly no solution.


restart:

L := [[1.0, 0.75e-1], [2.0, .1], [3.0, .12], [4.0, .14], [5.0, .16],
 [6.0, .18], [7.0, .2], [8.0, .22], [9.0, .24], [10.0, .26]]:

NN := -N+(N0B-NB)*(erf((1/2)*x/(sqrt(t)*sqrt(d)))
     -sqrt(erf((1/2)*alpha/d)))/sqrt(erfc((1/2)*alpha/d))+NB:

 

alpha:= 2*10^(-12):  NB:= 0.075:  N0B:= 0.2:  t:= 360000:

 

Digits:= 15:

for k to nops(L) do
     f:= unapply(eval(NN, [x,N]=~ L[k]), d);
     printf("\nTrying x= %f, N= %f\n", L[k][]);
     dd:= RootFinding:-NextZero(f, 0);
     if dd <> FAIL then
          printf("Solution found: d= %a\n", dd);
     else
          printf("No solution found. Showing plot instead.\n");
          print(plot(f, 0..1, labels= [d, 'NN']))
     end if
end do:


Trying x= 1.000000, N= 0.075000

No solution found. Showing plot instead.


Trying x= 2.000000, N= 0.100000
Solution found: d= .864547473017427e-4

Trying x= 3.000000, N= 0.120000
Solution found: d= .570969017295853e-4

Trying x= 4.000000, N= 0.140000
No solution found. Showing plot instead.


Trying x= 5.000000, N= 0.160000
Solution found: d= .350843030752835e-4

Trying x= 6.000000, N= 0.180000

Solution found: d= .253007839385318e-4

Trying x= 7.000000, N= 0.200000
No solution found. Showing plot instead.


Trying x= 8.000000, N= 0.220000
No solution found. Showing plot instead.


Trying x= 9.000000, N= 0.240000
No solution found. Showing plot instead.


Trying x= 10.000000, N= 0.260000
No solution found. Showing plot instead.

> 

 


Download erf_solve.mw

NonlinearFit...

Carl Love 28095
January 26 2014
0 9

Use Statistics:-NonlinearFit. Your data is the listlist L, and your parameter is d. Your model function is NN+N (so that it is in solved-for-N form).

You should change the title and tags of this Question. The fact that the data came from Excel, or that they were imported from anywhere, is totally irrelevant.

convert(..., hex)...

Carl Love 28095
January 24 2014
1 2

All that you need is

str1:= convert(Xor(f(number10), f(number8)), hex);

Note that there are no quotation marks in this command!

The command convert(..., bytes) has nothing to do with this. It deals with ASCII character values. However, if you still want to know the difference between the two forms of convert(..., bytes), it is this: They are inverses of each other. To convert the string expr to a list of ASCII decimal character values, use L:= convert(expr, bytes). To convert L back to the string expr, use convert(L, bytes). Note that quotation marks are only used to enclose literal string values: see ?symbol .

You can automate the conversion from/to hex so that it will work for any of the operators in Bits like this,

Hex:= Op-> convert(Bits[Op]((x-> convert(x, decimal, hex))~([_rest])[]), hex):

And then call it like this:

Hex(Xor, number8, number10);

 

Not currently possible...

Carl Love 28095
January 23 2014
1 0

It is not currently possible to access programmatically the directory name of the currently active worksheet.

sqrfree and Sqrfree...

Carl Love 28095
January 22 2014
2 8


f:= x^6-x^5-4*x^4+2*x^3+5*x^2-x-2:

g:= x^6-5*x^5+12*x^4-6*x^3-9*x^2+12*x-4:           

 

Over Q[x]:

sqrfree(f);

[1, [[x-2, 1], [x-1, 2], [x+1, 3]]]

sqrfree(g);

[1, [[x^6-5*x^5+12*x^4-6*x^3-9*x^2+12*x-4, 1]]]

Over Z3[x]:

Sqrfree(f) mod 3;

[1, [[x+2, 2], [x+1, 4]]]

Sqrfree(g) mod 3;

[1, [[x^6+x^5+2, 1]]]

 


Download Sqrfree.mw

help page...

Carl Love 28095
January 21 2014
0 0

There is some information pertaining to your question on the help page ?dsolve,numeric,BVP .

Iterator for loop...

Carl Love 28095
January 21 2014
2 8

I meant quitting the iterator programmatically. To do this, the iterator needs to be the index in a for loop.


restart:

You can make the second 1477 unique by making it 1477.00; there's no need to change the numeric value.

S:= [1829.0, 1644.0, 1594.0, 1576.0, 1520.0, 1477.0, 1477.00, 1404.0, 1392.0, 1325.0, 1313.0, 1297.0, 1292.0, 1277.0, 1249.0, 1236.0]:

SL:= [seq(A||i,i=1..nops(S))]:

assign(Labels ~ (S) =~ SL); #Create remember table.

 

lnp:= evalf(ln(`*`(S[])^(1/4))): #not ^(1/3)

Var:= proc(P::list(set))

local r:= evalf(`+`(map(b-> abs(ln(`*`(b[]))-lnp), P)[]));

end proc:

ts:= time():  
for P in Iterator:-SetPartitions({S[]},[[4,4]], compile= false)
     while Var(P) > .74 do
end do:

time()-ts;

 

subsindets(P, realcons, Labels);  

 

4.750

[{A14, A15, A16, A7}, {A13, A4, A5, A6}, {A11, A12, A2, A3}, {A1, A10, A8, A9}]

Var(P);

.6822751484438

 


Download Iterator_forloop.mw

sin(Pi/16)...

Carl Love 28095
January 21 2014
0 0


We get cos(Pi/16) by solving cos(8*x) = cos(Pi/2) = 0 for cos(x).

expand(cos(8*x))=0;

128*cos(x)^8-256*cos(x)^6+160*cos(x)^4-32*cos(x)^2+1 = 0

subs(cos(x)= C, %):

R:= [solve(%)];

[(1/2)*(-(2-2^(1/2))^(1/2)+2)^(1/2), -(1/2)*(-(2-2^(1/2))^(1/2)+2)^(1/2), (1/2)*(-(2+2^(1/2))^(1/2)+2)^(1/2), -(1/2)*(-(2+2^(1/2))^(1/2)+2)^(1/2), (1/2)*((2-2^(1/2))^(1/2)+2)^(1/2), -(1/2)*((2-2^(1/2))^(1/2)+2)^(1/2), (1/2)*((2+2^(1/2))^(1/2)+2)^(1/2), -(1/2)*((2+2^(1/2))^(1/2)+2)^(1/2)]

We take the largest root. To get sin(Pi/16), use sin(x) = sqrt(1-cos(x)^2), taking the positive branch.

sqrt(1-R[-2]^2);

(1/2)*(-(2+2^(1/2))^(1/2)+2)^(1/2)

evalf(% - sin(Pi/16));

0.2e-14

 


Download sinPi16.mw

Not thread safe...

Carl Love 28095
January 20 2014
1 2

LinearAlgebra:-Eigenvalues is not in the list at ?threadsafe, so there's no point trying to use it in multithreaded code.

Assignment is :=, not =...

Carl Love 28095
January 20 2014
2 0

Change rho[0]=1.33 to rho[0]:= 1.33. Likewise for T[0].

Applications Center...

Carl Love 28095
January 18 2014
2 1

A search for "quaternions" at the Maple Applications Center turns up a handful of packages.

with my logic package...

Carl Love 28095
January 18 2014
1 2

This type of problem can be done with my LogicProblem package, available at the Maple Applications Center, although in some ways it feels like using an elephant gun to swat a fly. One advantage over the other methods is that if the solution is unique, then my code will prove it unique.

 

restart:

read "C:/Users/Carl/desktop/logic_problems.mpl";

Vars:= ['Number','Letter']:

Number:= [$1..5]:  Letter:= [A,E,F,Z,P]:

LP:= LogicProblem(Vars):
with(LP);

Warning, LP is not a correctly formed package - option `package' is missing

[`&!!`, `&-`, `&<`, `&>`, `&?`, `&G`, `&Soln`, AutoGuess, CPV, CollectStats, ConstNum, Consts, ConstsInV, DifferentBlock, Equation, FreeGuess, GoBack, Guess, InternalRep, IsComplete, IsUnique, NC, NV, PrintConst, Quiet, Reinitialize, SameBlock, Satisfy, Separated, UniquenessProof, VarNum, VarNumC, X_O]

S:= [A > E, F > Z, F < P, Z > E, P < A]: #Input inequalities

Satisfy(subsindets(S,`<`, J-> Less(op(J), Number)));

NULL

`Unique solution:`

Matrix(5, 2, {(1, 1) = 1, (1, 2) = E, (2, 1) = 2, (2, 2) = Z, (3, 1) = 3, (3, 2) = F, (4, 1) = 4, (4, 2) = P, (5, 1) = 5, (5, 2) = A})

 

 

Download TopoSort.mw

First 323 324 325 326 327 328 329 Last Page 325 of 395