0
$\begingroup$

First step (done):

I start with my concrete problem: (sry for posting the full code, but might be better understanding)

Variables:

ρ0 := 1.2041;
P := 0.95;
G := 100000;
f := 6000;
ω := 2*Pi*f;
η := 1.85*10^-5;
j := 1.4;
P := 0.71;
B := Sqrt[P];
P0 := 101325;
(*k & z for f=6000*)
k = 170.09 - 82.22 I;
z = 676.6 - 330.61 I;

some simple equations:

α = (G*P)/(ρ0*ω)

k2 = Sqrt[z^2/(k/ω)^2]
ρ = z^2/k2

And the function with the searched variables "r" and "e":

Solve[ρ == ρ0*α*(1 + (P*G)/(
      I (ω*ρ0*α))*(1 + 
        I (4*ω*ρ0*η*α^2)/(G^2*P^2*r^2))^(
      1/2)), r]

Solve[k2 == P0/(
  j - (j - 1) (1 + (8*η)/(
       I (ω*B^2*e^2*ρ0))*(1 + 
         I (ω*B^2*ρ0*e^2)/(16*η))^(1/2))^-1), e]

My problem right now is that the variable f (=frequency) is continous (somewhere from 350 to 6400 in step 0.5).To every frequency there is a corresponding k and z.

I have my values for f , k and z in two .txt files. See at the end, there is a short ectract.

So, how to import these data correctly from the .txt-file (e.g. with space character between real and imaginary number) and solve, respectively, plot the searched variables "r" and "e" as a function of "f" with corresponding "z" and "k"?

I hope for concrete answers to my problem.

Thanks a lot!

f in Hz Re(k) in 1/m    IM(k) in 1/m

300 18,1021587151   -20,268568883
301 17,8270460036   -19,9416539488
302 17,574903357    -19,7650284049
303 18,8974277344   -19,344621652
304 18,799546006    -18,4230245471
305 19,3050841571   -18,0620751669
306 20,0889951843   -17,3210376494
307 21,1941412976   -17,0011098871

f in Hz Re(Z) in Ns/m^3 IM(Z) in Ns/m^3

300 1889,32930532   -1282,64805963
301 1924,55593429   -1283,96093068
302 1927,73749935   -1313,43266941
303 1846,3736042    -1331,9786319
304 1901,82478069   -1305,0076572
305 1903,96509609   -1310,07794373
306 1889,26900477   -1348,41152746
307 1853,51057624   -1344,98426415

Second step:

I want to inverse this method. So again my variables,

\[Rho]0 := 1.199;
\[Phi] := 0.97;  
\[Sigma] := 47742.675; 
\[Omega] := 2*Pi*f;
f = [(1000, 6300)];
\[Eta] := 1.81*10^-5;
\[Kappa] := 1.4;   
Pr := 0.72; 
B := Sqrt[Pr]; 
P0 := 98900;  
\[Lambda] = 1.296*10^-5;    
\[Xi] = 6.963*10^-5; 
\[Alpha] = 1;    
ko = \[Eta]/\[Sigma]; 
d = 14.1*10^-3;  
c0 = 345.53841;

But now, I want to solve following equations in dependency of frequency (f=1000-6300; how to define the steps?) AND for each frequency for different angle values [Epsilon] = [(0°, 90°)](define the single steps?). The aim is the average value of tau for each frequency of different angles.

Try of solve, export and plot the result.

result = Map[(f;
    z = Sqrt[k2*\[Rho]]; k = \[Omega]*Sqrt[\[Rho]/k2]; 
    Zs = z*coth[-I*k*d]; 
    R = (Zs/(\[Rho]0*c0)*cos (\[Epsilon]) - 1)/(Zs/(\[Rho]0*c0)*cos (\[Epsilon]) + 1); 
    \[Tau] = 1 - Abs[R]^2;
    {f, Solve[\[Rho] == \[Rho]0*\[Alpha]*(1 + (\[Phi]*\[Sigma])/(I (\
\[Omega]*\[Rho]0*\[Alpha]))*(1 + 
              I (4*\[Omega]*\[Rho]0*\[Eta]*\[Alpha]^2)/(\[Sigma]^2*\
\[Phi]^2*\[Lambda]^2))^(1/2)), \[Rho]], 
     Solve[k2 == 
       P0/(\[Kappa] - (\[Kappa] - 
             1) (1 + (8*\[Eta])/(I (\[Omega]*
                    B^2*\[Xi]^2*\[Rho]0))*(1 + 
                  I (\[Omega]*B^2*\[Rho]0*\[Xi]^2)/(16*\[Eta]))^(1/
                  2))^-1), k2], z, k, Zs, R, alpha})];

Export["result.txt", result];
Print[ListPlot[Map[{#[[8]]} &, result], 
   PlotLabel -> "\[Tau] versus f"]];

If there are again more solutions for the absolut value, there is only necessity of one of them.

Thanks a lot again!

$\endgroup$
2
$\begingroup$

I don't know whether in your location a comma is the default separator for decimal numbers. Here the period is the default separator. To make this code slightly simpler I have manually changed all the , to . in your two text files.

Do your variable initalizations and then

rawk = Import["kdata.txt", "Data"];
fmtk = Map[{#[[1]], #[[2]] + I*#[[3]]} &, rawk];
rawz = Import["zdata.txt", "Data"];
fmtz = Map[{#[[1]], #[[2]] + I*#[[3]]} &, rawz];
fmtkz = MapThread[Join, {fmtk, fmtz}];

reads the two data files, constructs each Complex from the two parts and finally joins the k and z for each row of your data file.

Then

result = Map[
   (f = #[[1]]; k = #[[2]]; z = #[[4]];
    α = (G*P)/(ρ0*ω); k2 = Sqrt[z^2/(k/ω)^2]; ρ = z^2/k2;
    {f, α, 
     Abs[r /. Solve[ρ == ρ0*α*(1 + (P*G)/(I (ω*ρ0*α))*(1 + 
       I (4*ω*ρ0*η*α^2)/(G^2*P^2*r^2))^(1/2)), r][[2]]], 
     Abs[e /. Solve[k2 == P0/(j - (j - 1) (1 + (8*η)/(I (ω*B^2*
       e^2*ρ0))*(1 + I (ω*B^2*ρ0*e^2)/(16*η))^(1/2))^-1), e][[2]]]
    })&, fmtkz];
Export["result.txt", result];
Print[ListPlot[{Map[{#[[1]], #[[3]]} &, result], 
   Map[{#[[1]], #[[4]]} &, result]}, PlotLabel->"Abs[r] and Abs[e] versus f"]];
Print[ListPlot[Map[{#[[1]], #[[2]]} &, result], PlotLabel -> "α versus f"]];

will take each row of your combined k and z, do the two Solve and give you the result.

EDIT: Revised to discard multiple roots, format and Export the results and ListPlot the results

$\endgroup$
  • $\begingroup$ Well, what should I say, that works perfectly! Thank you very much! $\endgroup$ – Sebastian Jun 7 '17 at 17:13
  • $\begingroup$ Well, what should I say, that works perfectly! Thank you very much! Theres is just this one point you already mentioned, one gets always two solutions for each k and z which are only different in the sign (+/-) My target is to plot the absolut value of "r" and "e" as well as the alpha-value as function of the frequency. So,it would be perfect, if I can 1) export the absolut value of "r" and "e" as well as alpha to a .txt file with corresponding frequency 2) plot the absolut values of "e" and "r" as well as alpha in mathematica If you could manage it too, that would be great! $\endgroup$ – Sebastian Jun 7 '17 at 17:27
  • $\begingroup$ @Sebastian Done. There is a lot going on in that code. You should go through that a character at a time studying that, taking it apart, trying to understand all the details of what has been done with that. Most of it is fairly simple, but there are some things in that which will not be obvious to a new MMA programmer. If you really understand every detail of how this works then you will have new skills to be better able to do your next task. And you should carefully test this to make certain there are no errors. The code is very dependent on your data being formatted exactly correctly. $\endgroup$ – Bill Jun 7 '17 at 19:00
  • $\begingroup$ Again Bill, I have to say that your work is excellent and very helpfully for my task of computing these variables. You did a great job and I could save a lot of time ! Thanks! The most things i could already understand, and I'm looking forward to learn more and more and getting better in mathematica! Thanks! $\endgroup$ – Sebastian Jun 8 '17 at 18:55

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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