void thermal_properties(void) { // // plot the number of fits to a cdc segment that result in prob > prob_cut (currently set to 0.01). // // #include #include gROOT->Reset(); //TTree *Bfield = (TTree *) gROOT->FindObject("Bfield"); gStyle->SetPalette(1,0); gStyle->SetOptStat(kFALSE); gStyle->SetOptFit(kFALSE); // gStyle->SetOptFit(1111); gStyle->SetPadRightMargin(0.15); gStyle->SetPadLeftMargin(0.15); gStyle->SetPadBottomMargin(0.15); gStyle->SetFillColor(0); // char string[256]; char filename[80]; Int_t j,jj; #define jpts 19; #define mpts 25; #define npts 3; #define max 201; # define npts_cp 10; Int_t ngen=10000; // temperature is in K, thermal conductivity in W/m-K, specific heat in J/kg-K // Double_t temperature[mpts]={1,2,3,4,5,6,7,8,9,10,15,20,30,40,50,60,70,80,90,100,150,200,250,300,350,400,500,600,800,1000,1200}; // Double_t conductivity[mpts]={4220,8400,12500,16200,19500,22200,23900,24800,24900,24300,17100,10800,4450,2170,1250,829,647,557,508,482,429,413,406,401,396,393,386,379,366,352,339}; Double_t temperature[mpts]={1,4,6,8,10,15,20,25,30,35,40,50,60,70,76,80,90,100,120,140,160,180,200,250,300}; Double_t conductivity[mpts]={-1,70,96,120,134,120,88,60,40,28,20,12,8,6.2,5.7,5.2,4.7,4.5,4.3,4.2,4.1,4.0,4.0,4.0,4.0}; Double_t conductivity_ofhc[mpts]={-1,2.4,3.7,4.7,6.0,8.5,11,12,12,11,10,7.7,6.2,5.5,5.2,4.9,4.7,4.5,4.3,4.2,4.1,4.0,4.0,4.0,4.0}; Double_t specific_heat[mpts]; Double_t thermal_diffusivity[mpts]; Double_t thermal_diffusivity_ofhc[mpts]; Double_t Cu_rho=8.96e3 ; // units are kg/m3 Double_t cu_cpt[npts_cp]={1,2,4,6,8,10,20,30,100,300}; Double_t cu_cp[npts_cp]={0.012,0.028,0.091,0.23,0.47,0.86,7.7,27,254,386}; // BNL Cryo p. VIII-B-2 J/kg-K Double_t etemp[jpts]={1,10,20,40,60,80,100,150,200,273,293,298,300,400,500,600,700,800,900}; Double_t electrical_resistivity[jpts]={0.002,0.00202,0.0028,0.0239,0.0971,0.215,0.348,0.699,1.046,1.543,1.678,1.712,1.725,2.402,3.09,3.792,4.514,5.262,6.041,}; Double_t L_WK = 2.44e-8; // Lorentz number, proportionality constant in the Wiedemann-Franz Law. Double_t electrical_resistivity_WK[mpts]; // calculation of resistivity from thermal condctivity and Wiedemann-Franz law. Double_t electrical_resistivity_WK_ofhc[mpts]; // calculation of resistivity from thermal condctivity and Wiedemann-Franz law. Double_t resistivity_formula[mpts]; // use formula from George's reference on MIITS calculation. (4.3-13) // TCanvas *c1 = new TCanvas("c1","c1 thermal_properties",200,10,700,700); c1->SetGridx(); c1->SetGridy(); c1->SetBorderMode(0); c1->SetFillColor(0); c1->SetLogx(); c1->SetLogy(); c1->Divide(2,2); c1->cd(1); c1_1->SetGridx(); c1_1->SetGridy(); c1_1->SetBorderMode(0); c1_1->SetFillColor(0); c1_1->SetLogx(); c1_1->SetLogy(); Double_t TD=343.5*0.85; // Debye temperature in Kelvin, emperical factor to reproduce data at about 4 K. Double_t xmin=0; Double_t xmax=500; Double_t ymin=0; Double_t ymax=10; // plot temperature function. Consider x as a parameter and plot time dependence. TF1 *debye = new TF1("debye_func",debye_func,xmin,xmax,1); Double_t x=0; // x is T/TD debye->SetParameter(0,TD); sprintf (string,"T_{D}=%.2f ^{o}K",TD); /* TF1 *debye1 = debye->DrawCopy("sameC"); leg->AddEntry(debye1,string,"l"); debye1->SetLineColor(2);*/ Double_t Cv; Double_t R=8.31; // Gas constant in J/Kmole Double_t Acu=63.546e-3; // Atomic weight of Cu in kg Double_t RRR=75; for (j=0;jIntegral(0,TDmax); // Cv = fudge*(R/Acu)*9*(T/TD)*(T/TD)*(T/TD)*sum; Cv = (R/Acu)*9*(T/TD)*(T/TD)*(T/TD)*sum; // printf ("j=%d T=%f TD=%f sum=%f, Cv=%f\n",j,T,TD,sum,Cv); specific_heat[j] = Cv; conductivity[j]=100*conductivity[j]; // convert from W/cm-K to W/m-K conductivity_ofhc[j]=100*conductivity_ofhc[j]; // convert from W/cm-K to W/m-K // conductivity_ofhc[j]=50*conductivity_ofhc[j]; // convert from W/cm-K to W/m-K and /2 for magneto resistive effect thermal_diffusivity[j] = conductivity[j]/(specific_heat[j]*Cu_rho); // use pure copper numbers thermal_diffusivity_ofhc[j] = conductivity_ofhc[j]/(specific_heat[j]*Cu_rho); // OFHC copper numbers electrical_resistivity_WK[j] = L_WK*T*1e8/conductivity[j]; // units of 10^-8 electrical_resistivity_WK_ofhc[j] = L_WK*T*1e8/conductivity_ofhc[j]; // units of 10^-8 resistivity_formula[j] = 1.545/RRR*(1+ 1/(2.32547e17*pow(T,-5) + 9.57137e13*pow(T,-3) + 1.62735e10*T)); printf (" T=%g, heat=%g, resistivity=%g, diffusivity=%g formula=%g\n",T,conductivity[j],electrical_resistivity_WK[j],thermal_diffusivity[j],resistivity_formula[j]); } // fudge on last point to adjust y scale // conductivity[mpts-1] = 100; Double_t ymin=100; Double_t ymax=20000; TGraph *conduct = new TGraph (mpts,temperature,conductivity); TGraph *conduct_ofhc = new TGraph (mpts,temperature,conductivity_ofhc); conduct->SetTitle(""); conduct->GetXaxis()->SetRangeUser(xmin,xmax); conduct->GetYaxis()->SetRangeUser(ymin,ymax); conduct->GetXaxis()->SetTitleSize(0.04); conduct->GetYaxis()->SetTitleSize(0.04); conduct->GetYaxis()->SetTitleOffset(1.5); // conduct->GetYaxis()->SetTitleOffset(2.0); conduct->GetXaxis()->SetTitle("Temperature (K)"); conduct->GetYaxis()->SetTitle("Thermal Conductivity k (W/m-K)"); conduct->GetXaxis()->SetNdivisions(505); conduct->SetMarkerColor(1); conduct->SetMarkerStyle(20); conduct->SetMarkerSize(0.8); conduct->Draw("Ap"); conduct_ofhc->SetMarkerColor(2); conduct_ofhc->SetMarkerStyle(20); conduct_ofhc->SetMarkerSize(0.8); conduct_ofhc->Draw("Samep"); // plot k from empirical resistivity and the Wiedemann-Franz law Double_t thermal_empirical[1]; Double_t temp_empirical[1]={4}; Double_t T = temp_empirical[0]; Double_t res = 1.7*(2/75.); // units if 10-8, room temperature value x magneto resistive effect / RRR. thermal_empirical[0] = L_WK*T*1e8/res; TGraph *thermal_emp = new TGraph(1,temp_empirical,thermal_empirical); thermal_emp->SetMarkerColor(2); thermal_emp->SetMarkerStyle(21); thermal_emp->Draw("Samep"); printf ("T=%g, L_WK=%g, res=%g, thermal_emp=%g\n",T,L_WK,res,thermal_empirical[0]); // // TCanvas *c2 = new TCanvas("c2","c2 thermal_properties",200,10,700,700); c1->cd(2); c1_2->SetGridx(); c1_2->SetGridy(); c1_2->SetBorderMode(0); c1_2->SetFillColor(0); c1_2->SetLogx(); c1_2->SetLogy(); TLegend *legcp = new TLegend(0.30,0.15,0.85,0.25); TGraph *specific = new TGraph (mpts,temperature,specific_heat); specific->SetMarkerColor(2); specific->SetMarkerStyle(20); specific->SetMarkerSize(0.8); specific->SetTitle(""); specific->GetXaxis()->SetRangeUser(xmin,xmax); // specific->GetYaxis()->SetRangeUser(ymin,ymax); specific->GetXaxis()->SetTitleSize(0.04); specific->GetYaxis()->SetTitleSize(0.04); specific->GetYaxis()->SetTitleOffset(1.5); specific->GetXaxis()->SetTitle("Temperature (K)"); specific->GetYaxis()->SetTitle("Specific Heat c_{p} (J/kg-K)"); specific->GetXaxis()->SetNdivisions(505); legcp->AddEntry(specific,"Debye Formula Cu","p"); specific->Draw("Ap"); TGraph *bnl_cp = new TGraph(npts_cp,cu_cpt,cu_cp); bnl_cp->SetMarkerColor(4); bnl_cp->SetMarkerStyle(20); bnl_cp->SetMarkerSize(0.8); bnl_cp->Draw("Samep"); legcp->AddEntry(bnl_cp,"BNL Data Cu","p"); legcp->Draw(""); // // TCanvas *c3 = new TCanvas("c3","c3 thermal_properties",200,10,700,700); c1->cd(3); c1_3->SetGridx(); c1_3->SetGridy(); c1_3->SetBorderMode(0); c1_3->SetFillColor(0); c1_3->SetLogx(); c1_3->SetLogy(); Double_t ymin=0.0001; Double_t ymax=20; TGraph *diffusivity = new TGraph (mpts,temperature,thermal_diffusivity); TGraph *diffusivity_ofhc = new TGraph (mpts,temperature,thermal_diffusivity_ofhc); diffusivity->SetTitle(""); diffusivity->GetXaxis()->SetRangeUser(xmin,xmax); diffusivity->GetYaxis()->SetRangeUser(ymin,ymax); diffusivity->GetXaxis()->SetTitleSize(0.04); diffusivity->GetYaxis()->SetTitleSize(0.04); diffusivity->GetYaxis()->SetTitleOffset(1.5); diffusivity->GetXaxis()->SetTitle("Temperature (K)"); diffusivity->GetYaxis()->SetTitle("Thermal Diffusivity #alpha (m^{2}/s)"); diffusivity->GetXaxis()->SetNdivisions(505); diffusivity->SetMarkerColor(1); diffusivity->SetMarkerStyle(20); diffusivity->SetMarkerSize(0.8); diffusivity->Draw("Ap"); diffusivity_ofhc->SetMarkerColor(2); diffusivity_ofhc->SetMarkerStyle(20); diffusivity_ofhc->SetMarkerSize(0.8); diffusivity_ofhc->Draw("Samep"); // // TCanvas *c4 = new TCanvas("c4","c4 thermal_properties",200,10,700,700); c1->cd(4); c1_4->SetGridx(); c1_4->SetGridy(); c1_4->SetBorderMode(0); c1_4->SetFillColor(0); c1_4->SetLogx(); c1_4->SetLogy(); TLegend *leg = new TLegend(0.15,0.85,0.85,1.0); TGraph *resistivity = new TGraph (jpts,etemp,electrical_resistivity); TGraph *resistivity_WK = new TGraph (mpts,temperature,electrical_resistivity_WK); TGraph *resistivity_WK_ofhc = new TGraph (mpts,temperature,electrical_resistivity_WK_ofhc); Double_t resistivity_empirical[1]={1.7*2/75}; // units if 10-8, room temperature value x magneto resistive effect / RRR. Double_t temp_empirical[1]={4}; TGraph *resistivity_emp = new TGraph(1,temp_empirical,resistivity_empirical); TGraph *resistivity_for = new TGraph(mpts,temperature,resistivity_formula); Double_t ymin=0.001; Double_t ymax=5; resistivity_WK->SetTitle(""); resistivity_WK->GetXaxis()->SetRangeUser(xmin,xmax); resistivity_WK->GetYaxis()->SetRangeUser(ymin,ymax); resistivity_WK->GetXaxis()->SetTitleSize(0.04); resistivity_WK->GetYaxis()->SetTitleSize(0.04); resistivity_WK->GetYaxis()->SetTitleOffset(1.5); resistivity_WK->GetXaxis()->SetTitle("Temperature (K)"); resistivity_WK->GetYaxis()->SetTitle("Electrical Resistivity_WK #rho (10^{-8}#Omega-m)"); resistivity_WK->GetXaxis()->SetNdivisions(505); resistivity_WK->SetMarkerColor(4); resistivity_WK ->SetMarkerStyle(20); resistivity_WK->SetMarkerSize(0.8); resistivity_WK->Draw("Ap"); sprintf (string,"k & Wiedemann-Franz Law Pure Cu\n"); leg->AddEntry(resistivity_WK,string,"p"); resistivity_WK_ofhc->SetMarkerColor(2); resistivity_WK_ofhc ->SetMarkerStyle(20); resistivity_WK_ofhc->SetMarkerSize(0.8); resistivity_WK_ofhc->Draw("Samep"); sprintf (string,"k & Wiedemann-Franz Law OFHC\n"); leg->AddEntry(resistivity_WK_ofhc,string,"p"); /*resistivity_for->SetMarkerColor(1); resistivity_for->SetMarkerStyle(20); resistivity_for->Draw("Samep"); sprintf (string,"Formula 4.3-13\n"); leg->AddEntry(resistivity_for,string,"p");*/ resistivity->SetMarkerColor(1); resistivity->SetMarkerStyle(20); resistivity->SetMarkerSize(0.8); resistivity->Draw("Samep"); sprintf (string,"Data on pure Cu\n"); leg->AddEntry(resistivity,string,"p"); resistivity_emp->SetMarkerColor(2); resistivity_emp->SetMarkerStyle(21); // resistivity_emp->SetMarkerSize(0.8); resistivity_emp->Draw("Samep"); sprintf (string,"Scaled from Room T by 2/75\n"); leg->AddEntry(resistivity_emp,string,"p"); leg->Draw(); Double_t alpha=1; sprintf(string,"#alpha/k =%.3f Kcm3/J\n",alpha); t1 = new TLatex(0.50,0.65,string); t1->SetTextColor(1); t1->SetNDC(); t1->SetTextSize(0.025); // t1->Draw(); sprintf(filename,"thermal_properties_c1.eps"); c1->SaveAs(filename); sprintf(filename,"thermal_properties_c1.png"); c1->SaveAs(filename); /*sprintf(filename,"thermal_properties_c2.pdf"); c2->SaveAs(filename); sprintf(filename,"thermal_properties_c2.png"); c2->SaveAs(filename); sprintf(filename,"thermal_properties_c3.pdf"); c3->SaveAs(filename); sprintf(filename,"thermal_properties_c3.png"); c3->SaveAs(filename); sprintf(filename,"thermal_properties_c4.pdf"); c4->SaveAs(filename); sprintf(filename,"thermal_properties_c4.png"); c4->SaveAs(filename);*/ } Double_t debye_func (Double_t *x, Double_t *par) { // Compute function inside integral in Debye function Double_t TD=par[0]; Double_t x1=x[0]; Double_t pi=3.14159; char string[256]; Double_t func; if (x1 < 0.0001) { func = x1*x1; } else if { func = (x1*x1*x1*x1*exp(x1))/((exp(x1)-1)*(exp(x1)-1)); } /*sprintf (string,"x1=%f func=%f\n",x1,func); printf ("string=%s",string);*/ return func; }