Categories
Mathematica

The dynamics of Plasmodium falciparum during the expansion phase of the asexual stage of infection

The dynamics of Plasmodium falciparum during the expansion phase of the asexual stage of infection

https://www.wolframcloud.com/objects/sompob/models/NJWandImmune.nb

 

Categories
Mathematica

เขียนcodeแก้ปัญหาจากpantipเล่น

https://pantip.com/topic/38221485/

code สำหรับปัญหานี้

https://www.wolframcloud.com/objects/sompob/test/pantip_01Nov18.nb

Manipulate[
Show[
  Graphics[{EdgeForm[Pink], FaceForm[Blue], 
    Rectangle[{1 - Cos[-15 Degree], 0}, {x, 
      x - (1 - Cos[-15 Degree])}]}],
  Graphics[{EdgeForm[LightBlue], FaceForm[Black], 
    Rotate[Rectangle[{0, 0}], -15 Degree, {1, 0}],
    PlotRange -> {{-2, 2}, {-2, 2}} }],
  Graphics[
   Text["Area of the Big Square:" <> 
     ToString[(x - (1. - Cos[15 Degree]))^2/(0.5*Sin[15 Degree]*
         Cos[15 Degree])], {0.25, 1.5}]]
   ],
 {{x, 1.259, "side"}, 0.1, 2, 0.001}
 ]

 

 

Categories
Mathematica

My Green Apple

ทดลองเล่น WolframCloud

เจ๋งมากกกกกก เขียนสั้น deploy ก็ง่าย (เชียร์ให้คนมาใช้กันเยอะๆ 555)

 

Categories
Mathematica R

แสดง output จาก ggplot ใน Mathematica

ลองดูตัวอย่างที่ผมเขียนนี้ครับ  ไอเดียก็คือว่าให้มันเขียนภาพเป็นไฟล์ pdf ด้วยคำสั่ง ggsave แล้วเรียกกลับมาใน Mathematica ด้วย Import  โค้ดที่เขียนนี้ดัดแปลงมาจาก https://mathematica.stackexchange.com/questions/16726/how-to-get-a-plot-generated-by-r-returned-in-an-output-cell

(* using R that already installed on the computer *)
R363 = "C:\\Data-Work\\Programs\\R-3.6.3\\App\\R-Portable";
InstallR["RHomeLocation" -> R363, "RVersion" -> "3.6.3", 
 "NativeLibLocation" -> 
  "C:\\Data-Work\\Programs\\R-3.6.3\\App\\R-Portable\\library\\rJava\\\
jri\\x64"]

REvaluate["
 library(ggplot2)
 "]
(* export ggplot as pdf file *)
Wrapper = RFunction["function(filename,plotfun){
   ggsave(filename, plot=plotfun(), width=6, height=4, \
device=cairo_pdf)
   }"];

(* generate the plot as pdf and import back to show in Mathematica *)


Getggplot[plotFun_] := Module[{tempfile, rcode, runcode},
   tempfile = FileNameJoin[{$TemporaryDirectory, "temp.pdf"}];
   If[FileExistsQ[tempfile], Quiet@DeleteFile[tempfile]];
   rcode = "ADA <- ggplotcmd;";
   runcode = StringReplace[rcode, {"ggplotcmd" -> plotFun}];
   Wrapper[tempfile, RFunction["function()" <> "{" <> runcode <> "}"]];
   If[! FileExistsQ[tempfile], Return[$Failed]];
   Import[tempfile] // First
   ];
Categories
Mathematica

มาลองทำให้โค้ดMathematicaรันเร็วขึ้น

โค้ดที่เขียนด้วยMathematicaหรือภาษา Wolfram นั้นมันจะมีวิธีที่ช่วยให้มันรันได้เร็วขึ้นอยู่หลายวิธีครับ ขึ้นกับปัญหาและวิธีการเขียน ที่จะโชว์ให้ดูนี้ผมก็ใช้เทคนิคง่ายๆด้วยการทำ parallel ด้วยคำสั่ง ParallelTable และ compile ด้วยคำสั่ง Compile ครับ  โค้ดตัวอย่างที่จะเอามาลองนี้เป็นโค้ดที่ใช้วาดรูป Mandelbrot set ครับ

แบบที่ 1 เป็นโค้ดเริ่มต้นที่ยังไม่มีการทำให้มันเร็วครับเขียนโดยใช้ loop จากคำสั่ง Table

AbsoluteTiming[
 test1 = Block[{i, x, p}, 
    Table[i = 0; x = 0. I; p = r + I c; 
     While[Abs@x <= Sqrt[2] && i < 9^3, x = x^2 + p; ++i]; 
     Tanh[Power[i/9^3, (7)^-1]], {c, -1, 1, .01}, {r, -2, 1, .01}]];]

(* เวลาที่ใช้คือ 77.4437 วินาที *)

แบบที่ 2 ในแบบนี้ผมทดลองให้มันรันแบบขนานโดยใช้คำสั่ง ParallelTable เครื่องที่ผมใช้นี้มันมี 4 cores ครับ เวลาก็ลดลง ถ้าคิดง่ายๆก็ใช้เวลาเหลือเกือบๆ 1 ใน 4

AbsoluteTiming[
 test2 = Block[{i, x, p}, 
    ParallelTable[i = 0; x = 0. I; p = r + I c; 
     While[Abs@x <= Sqrt[2] && i < 9^3, x = x^2 + p; ++i]; 
     Tanh[Power[i/9^3, (7)^-1]], {c, -1, 1, .01}, {r, -2, 1, .01}]];]

(* เวลาที่ใช้คือ 23.8938 วินาที *)

แบบที่ 3  แบบนี้ก็คือเอาแบบที่ 1 มา compile เป็น binary code ด้วยคำสั่ง Compile เวลาก็เร็วขึ้นอย่างเห็นได้ชัด

test3 = Compile[{}, 
  Block[{i, x, p}, 
   Table[i = 0; x = 0. I; p = r + I c; 
    While[Abs@x <= Sqrt[2] && i < 9^3, x = x^2 + p; ++i]; 
    Tanh[Power[i/9^3, (7)^-1]], {c, -1, 1, .01}, {r, -2, 1, .01}]]]

AbsoluteTiming[test3[];]

(* เวลาที่ใช้คือ 2.68316 วินาที *)

แบบที่ 4  แบบนี้ก็คือแบบที่ 3 ที่เอามารันแบบขนานด้วยคำสั่ง ParallelTable ซึ่งเวลาที่ใช้ก็เร็วขึ้นอีกกว่าแบบที่ 3

test4 = Compile[{{c, _Real}, {r, _Real}},
   Module[{i, x, p},
    i = 0; x = 0. I; p = r + I c; 
    While[Abs@x <= Sqrt[2] && i < 9^3, x = x^2 + p; ++i]; 
    Tanh[Power[i/9^3, (7)^-1]]
    ]
   ];

ParallelTable[
   test4[c, r], {c, -1, 1, .01}, {r, -2, 1, .01}]; // AbsoluteTiming

(* เวลาที่ใช้คือ 1.04455 วินาที *)

แบบที่ 5 แบบนี้เหมือนกับแบบที่ 4 แต่แทนที่จะให้มัน compile เป็น binary code โดย Mathematica เอง ก็ให้มันแปลงเป็นภาษา C ด้วย option CompilationTarget -> “C” เลยโดยใช้ compiler อย่าง gcc ครับ เวลาก็เร็วขึ้นอีกประมาณ 60% จากแบบที่ 4 ครับ

Needs["CCompilerDriver`GenericCCompiler`"]

$CCompiler = {"Compiler" -> GenericCCompiler, 
   "CompilerInstallation" -> "C:/Rtools/mingw_64", 
   "CompilerName" -> "x86_64-w64-mingw32-gcc.exe"};


test5 = Compile[{{c, _Real}, {r, _Real}},
   Module[{i, x, p},
    i = 0; x = 0. I; p = r + I c; 
    While[Abs@x <= Sqrt[2] && i < 9^3, x = x^2 + p; ++i]; 
    Tanh[Power[i/9^3, (7)^-1]]
    ], CompilationTarget -> "C"];

ParallelTable[
   test5[c, r], {c, -1, 1, .01}, {r, -2, 1, .01}]; // AbsoluteTiming

(* เวลาที่ใช้ 0.393584 วินาที *)

 

Categories
Mathematica

external packages / tools / resources ของ Mathematica ที่สนใจ

แนะนำ Link นี้ครับสำหรับผู้สนใจใช้งาน Mathematica ขุมทรัพย์อีกแห่งเลยครับ 🙂

https://stackoverflow.com/questions/4198961/what-is-in-your-mathematica-tool-bag

เทคนิคเจ๋งๆเพียบเลยครับ

 

 

Categories
Mathematica R

ใช้งาน gcc ที่มากับ Rtools ใน Mathematica

ใน Rtools จะมี compiler ของ gcc (mingw) มาด้วยแล้วทั้งที่เป็นแบบ 32 bits และ 64 bits หากอยากจะเอาไปใช้ใน Mathematica ก็ต้องเรียกผ่าน CCompilerDriver`GenericCCompiler` ครับ แล้วเพียงเซ็ต Path ของ gcc จาก Rtools นี้ให้ถูก เช่น

Needs["CCompilerDriver`GenericCCompiler`"]

ทดลองเรียกใช้งาน

greeter = CreateExecutable[StringJoin[
 "#include <stdio.h>\n",
 "int main(){\n",
 " printf(\"Hello MinGW-w64 world.\\n\");\n",
 "}\n"],
 "helloworld", "Compiler" -> GenericCCompiler, 
 "CompilerInstallation" -> "C:/Rtools/mingw_64", 
 "CompilerName" -> "x86_64-w64-mingw32-gcc.exe"]

Import["!\""<>greeter<>"\"","Text"]
Hello MinGW-w64 world.

ที่นี้ถ้าอยากจะให้ Mathematica มันเรียกใช้เจ้า mingw 64 จาก Rtools นี้ตลอด อย่างเช่นคำสั่ง Compile ก็สามารถทำได้โดยเซ็ตค่าที่มันเกี่ยวข้องอย่างเช่น Path กับตัวแปร $CCompiler ได้เลยครับ

$CCompiler = {"Compiler" -> GenericCCompiler, 
 "CompilerInstallation" -> "C:/Rtools/mingw_64", 
 "CompilerName" -> "x86_64-w64-mingw32-gcc.exe"};

f = Compile[{x, y}, Sqrt[x^2 + y^2], CompilationTarget -> "C"]

Table[{x, f[x, 5/4 x]}, {x, 0, 6, 0.5}]

 

Categories
Mathematica R Uncategorized

ทดลองใช้ C# กับงานวิจัย

ปกติผมจะใช้ Mathematica เป็นหลักในการทำงานโดยมากจะเป็นพวกแบบจำลองคณิตศาสตร์ที่ตองการเห็นว่าผลลัพธ์จะเปลี่ยนแปลงไปอย่างไรถ้าตัวแปรที่สนใจบางตัวมีการเปลี่ยนแปลง เพิ่มขึ้นหรือลดลง ซึ่งแน่นอนผมใช้พวกคำสั่ง Manipulate หรือ Dynamic เป็นหลัก แต่มาระยะหลังนี้เริ่มมีการใช้ R ในที่ทำงานกันเยอะ ผมเลยเขียนโมเดลของงานที่ทำใน R ไปพอสมควร แต่มันก็ติดปัญหาว่ามันช้ามากกับหลายโมเดลที่ทำอยู่ถึงแม้จะหันไปใช้ Rcpp แต่ก็รู้สึกว่าเสียเวลาเขียน wrapper functions ห่อพวกมันอีกทีใน R แถมเวลาส่งงานให้ Prof ทั้งหลายก็ติดเรื่อง compiler อีก แถมต้องมาเขียน shiny อีกดูวุ่นวายไปหมด สุดท้ายกลับมาตายรังที่ Mathematica เหมือนเดิม เขียนโมเดลเสร็จ อยากจะส่งให้ชาวบ้านดูก็ใช้เพียงแค่คำสั่ง Deploy ไปบนcloud เท่านั้นจบเลย …แต่มันก็แลกมาด้วยค่าใช้จ่ายพอสมควร

ผมเลยหาทางเลือกอื่นๆที่มันประหยัดกับบางงานกับโมเดลที่ไม่ดูซับซ้อนและใช้เวลาทำไม่นาน หันซ้ายหันขวาก็มาจบที่ C# หลังจากดูแล้วก็มี library ที่พอช่วยงานได้อย่าง math.net แถม UI ก็ทำง่ายเพียงลากแปะเท่านั้น เลยทดลองใช้ดู ผลปรากฏว่าเป็นที่น่าพอใจ  ถึงจะดูยุ่งยากในส่วน UI แต่ก็รู้สึกว่าไม่ได้เลวร้ายอะไร ส่งงานก็ง่ายแถมก็ไม่ต้องติดตั้งอะไรให้ดูวุ่นวาย

ตัวอย่างงงานที่ทำที่พอให้ดูได้ครับ


สนใจอยากทดลองเล่นก็ download ได้เลยครับ

Categories
Mathematica

ปัญหา cuTexRefSetAddress2D_v3

พอดีว่าลงWindows 7 64bit กับ Mathematica 11.1.1 ใหม่ แล้วเจอปัญหา  ตามภาพด้านล่างตอนใช้ CUDALink พอเช็คดูเวอร์ชั่นของ nvidia driver ที่ใช้กับการ์ดจอ (Quadro 2000)ด้วยคำสั่ง CUDADriverVersion[] พบว่าเป็นรุ่นเก่ามากคือ27x.x  เลยทำการ update เป็นตัวใหม่กว่าคือ 377.55จากนั้นก็ลองใช้ QUDAQ[] ใหม่เช็คดูปรากฎว่าไม่มีปัญหาแล้วสามาถใช้งานได้ตามปกติ

Categories
Mathematica

ปัญหาเข็มของบึฟฟองBuffon’s Needle

เห็นมีคนพูดถึงปัญหานี้ในกลุ่มผู้ใช้ R ของไทยเราเลยเอามาเขียนโปรแกรมทำบ้างแต่ใช้ภาษา Wolfram แทนเพราะจะเขียนได้สั้นมาก

โปรแกรมที่เขียนนี้จะวาดเส้นแนวนอนที่มีระยะห่างตามที่กำหนดและก็สุ่มวางเข็มตามขนาดและจำนวนที่ต้องการลงไปแล้วนับดูว่าที่กี่เล่มที่ทับเส้น จากนั้นก็คำนวณประมาณค่า pi จากสูตรที่เขียนไว้ที่ https://en.wikipedia.org/wiki/Buffon%27s_needle ซึ่งก็คือ

โดย l คือความยาวเข็ม, n  คือจำนวนเข็ม, t คือระยะห่างระหว่างเส้นบรรทัด และ h คือจำนวนเข็มที่ทับเส้นบรรทัด

สำหรับวาดเส้นของเข็มและเส้นบรรทัด


genNeedle[{x_, y_}, θ_, len_] :=
Line[{{(x + 0.5 len Cos[θ]), (y + 0.5 Sin[θ])}, {(x –
0.5 len Cos[θ]), (y – 0.5 len Sin[θ])}}]

horizontalLine[diff_, len_, {x0_, y0_}, n_] :=
Table[Line[{{x0, y0 + (i diff)}, {len, y0 + (i diff)}}], {i, n}]

สำหรับประมาณค่า pi


estimatedPi[needlelength_, nneedles_, linespace_, ncross_] :=
2.*needlelength*nneedles /(linespace*ncross)

สำหรับ simulation


run[nl_, nfl_, ll_, diff_] :=
Module[{fls, ls, outls, pnts, ncross},
fls = horizontalLine[diff, 10, {0, 0}, nfl];
ls = Table[
genNeedle[{RandomReal[{0, 10}], RandomReal[{0, nfl*diff}]},
RandomReal[π], ll], {nl}];
outls =
Table[Solve[{x, y} ∈ RegionIntersection[fls[[i]], #], {x,
y}] & /@ ls, {i, Length@fls}];
pnts = Point[{x, y}] /. Flatten[outls, 2];
ncross = Length@Flatten[outls, 2];
Graphics[{fls, ls, Red, PointSize[Medium],
If[Length@Flatten[outls, 2] > 0, pnts]},
PlotLabel ->
Style["จำนวนเข็มที่ทับเส้น:" <> ToString@ncross <>
" π≈" <>
ToString@estimatedPi[ll, nl, diff, ncross], FontSize -> 20,
Bold]]];

ทดลองเล่นดู โดยใช้เข็ม 100 เล่ม และกำหนดให้มีเส้นบรรทัด 10 เส้น โดยระยะห่างระหว่างบรรทัดกับความยาวเข็มเท่ากัน


nlines = 10;
nneedles = 100;
linelength = 1;
linespace = 1;
run[nneedles, nlines, linelength, linespace]

สำหรับmathematica notebook สามารถdownloadได้จากhttps://github.com/slphyx/Buffon-s-Needle