Photocatalysis

Photocatalysis is the acceleration of photoreaction in the presence of a catalysis, with the help of light.  It can find many applications such as

• Clean energy: split water to hydrogen
• Sterilization: sterilize water or surgical instruments
• Distillation & desalination: desalinate seawater through solar vapor generation
• Decomposition: Photodegradate environmentally harmful matters  such as dyes, organic compounds,
• Photosynthesis: synthesize useful chemicals from environmentally harmful chemicals such as converting CO2 with H2O to CH4
• bio-medical: Photomelt, photoweld, halogenate DNA nucleobases

This page introduces what we can and cannot simulate as well as simulation methodology. Detailed simulation examples and analysis are also provided.

Introduction

To enhance the photocatalysis process, recently the excited localized surface plasmonics (LSP) in metallic nanostructures attracts much attention, and some authors term it as thermal plasmonic photocatalysis. With proper design, the illuminating light power can be coupled to the inhomogeneous oscillating electric fields in the nanostructure, and the induced intensified electric fields can increase the light absorption due to dissipation. From the energy point of view, the resonant state at high energy can decay in the forms of re-emission (scattering) or charge carrier. Mainly due to Joule effect, the absorbed light causes temperature to increase. Under CW illumination, the metallic nanoparticles (NPs) will thermally equilibrate internally quickly. The temperature raise can be spatially localized around the nano particles or smoothed due to collective thermal effects depending on such as the density of the NPs. In addition to plasmonic photocatalysis, semiconductor photocatalysis is also promising, and their combination may lead to new opportunities.

Photocatalysis involves optical, thermal, chemical and sometimes electrical effects. We are interested only in optical and thermal simulations and will not explore its physics.

Optical simulations: Such simulations of photocatalysis using FDTD are aimed to design the nanostructures in order to increase the capability of capturing incident light by tackling the properties of the nanostructures, such as the shape, size, orientation, inter-particle distance or arrangement, NP materials, surrounding materials etc to lead strong plasmonic resonance at the desired wavelengths or high absorption rate at a given wavelength range.

Thermal simulation: Such simulations of photocatalysis using HEAT Thermal Transport Solver is to find the temperature increase and thermal distribution in the nanostructures due to light absorption.

What we can simulate:

For optical simulations, the following results can be obtained (but not limited to):

total power absorbed (W), absorption density as a function of space and/or wavelength, extinction cross section or efficiency, absorption cross section, scattering cross section, reflection, transmission/reflecton, electrical field intensity profile, field enhancement and overall enhancement, near field confinement, plasmonic effect analysis such as dipole or quadrupole resonances etc. As the original outputs of FDTD Simulations are electric and magnetic fields, any quantities related to them can be obtained through post processing, as long as the physics is included in the Maxwell’ Equations.

For thermal simulations, the following results can be obtained (but not limited to):

Volumetric heat source density, or optical absorption density (W/m3), thermal gradient (K/m or K/um), thermal time constant (K/s), Thermal cross section, temperature change of spatial distribution, Temperature vs light flux.

For both optical and thermal simulations, transient and equilibrium state or steady state can be obtained.

What we do not simulate:

Charge transfer or photo-current except for semiconductors, vapor, photocatalytic enhancement, reaction rate, photocatalytic efficiency, chemical activities etc that are not related to Maxwell’s Equations and thermal transfer equations.

Simulation methodology

The photo-thermal-electrical-chemical are occurring concurrently, but we simulate them separately. Thus some material properties may change due to the complicated processes, such as refractive index and specific heat coefficient. The refractive index change due to temperature increase usually is small for most materials since the thermal coefficient is generally in the range of 10-4 ~10-7/K. Therefore the optical simulation results may not usually be sensitive to the temperature in photcatalysis. However, optical material properties may need to be modified due to some quantum effects, eg when two nano metallic particles with the same material are close enough that their permittivity may need to be modified from their bulk property to take the quantum effect into consideration.

Material’s thermal property can be sensitive to the geometric dimension. The bulk thermal conductivity is the upper limit. On the nano scale, the thermal conductivity can be reduced to as low as 10% of their bulk value. So care must be taken to modify them based on experiment or literature. may depend on temperature and can be anisotropic.

Material’s thermal property such as thermal conductivity can also be sensitive to the temperature, which can be turned on in the material thermal property.

Furthermore, users need to consider how the heat is transferred from NPs to environment. For example, in distillation and desalination, the heat produced from the NPs is transferred directly to water and vaporize it in very small region with boiling temperature, and some amount of heat can be conducted into substrate. However the result could be wrong if the thermal impedance of the substrate or the solver boundary condition is not set properly.

Optical simulations

There are four fundamental steps in setting up the simulation project file:

Physical structure

geometric shapes can be fundamental forms such as sphere, shell, ellipsoid, triangle/ pyramid/cube/ rod /disk/ring /slab or combinations such as bowtie, dimer, trimer etc. The arrangement can be single, arrayed, or finitely arrayed. Interested materials are some noble metals, semiconductors, dielectrics and some alloys, and the mostly used materials can be found in FDTD’ Material Library. If not there, users can import material data.

• Size: The geometrical size of NPs is typically from 1 to 100nm but can be as large as 10 microns, depending on purposes.
• Single or coupled array: in practice single nano particle is rarely used, however it may be helpful to understand the mechanisms and some properties such as extinction efficiency, enhancement of plasmonic resonance etc. It has been shown that gaped nano particles can have orders of magnitudes of electric field enhancement. Thus arrayed NPs are commonly used. In addition, it significantly increases the volume of interest. Square, hexagonal or cloud-particle array are usually used, and in some cases the NPs are densely packed. In some experiments laser beam is focused on a small region of array such as 3 to 4 periods with other considerations such as easy vaporization and flow. The structure groups in the Object Library can be used for complicated structures.
• Randomness: Chemically-grown nanoparticles may be difficult to control size/shape/pitch, and nanoparticles may be randomly aggregated on a surface; Other methods will inevitably have manufacturing errors. Therefore different sizes and spatial distributions may need to be simulated statistically. The object library in FDTD has some structure groups that can be used or modified. Such statistical simulations can be very challenge as many statistical parameters such as correlation length, probability function (Gaussian, uniform, triangle etc), mean values and standard deviation as well as density can be difficult to obtain accurately.

For each geometrical object, its material optical property of refractive index or permittivity must be assigned from Material Database. Before simulation, users should check or/and adjust the material fittings in the given spectra bandwidth.

FDTD uses mesh order to prioritize overlapped geometries, which can be modified in the Material Tab of the object.

Simulation region, Boundary conditions and Mesh

Simulation regions are limited by proper boundary conditions. Boundary conditions may depend on source type and geometry. Finer override mesh can be used for efficient simulation. Although the transmission and reflection results can be improved with finer mesh and/or higher mesh accuracy, the resulting temperature may not improve much.Users are encouraged to test using proper yet efficient mesh size and accuracy.

Optical source

The most often used source for photocatalysis simulation is plane wave and TFSF. Sometimes a focused Gaussian beam can also be used. For normal incidence, plane waves can be used to simulate array together with periodic BCs. TFSF can be used to simulate single NPs such as sphere or dimmer together with PML BCs. Sometimes TFSF can also be used for periodic structures.

• Source spectrum: In experiment, solar light (broadband) or laser light (single frequency) can be used. In simulations, it is not necessary to use real light source. Instead the built-in broadband source can be used as long as the users set the correct wavelength range. This is because in photocatalysis it usually involves linear optical interactions and the result is normalized to source spectrum. The real light absorption is obtained through post processing by multiplying the normalized result with actual light spectrum such as solar.
• Source power: The absorption density calculated in the analysis group is proportional to the input source power. So in optical simulation usually the incident power is not important. However the real source power can be scaled in thermal simulation when the absorption density is imported to thermal solver. Due to normalization the default source amplitude of 1 V/m can be used.
• Pulsed or CW: Pulsed source is always used in linear optics. In most cases the CW result is desired for optical simulations. The time scale for photocatalysis is in the order of milli- or nano- seconds however in simulation it is on the order of nano- or femto- seconds. The pulse source is most widely used yet the CW results can be obtained through frequency-domain monitors, in the linear optics regime.
• Source polarization: SPR is sensitive to polarization. For solar light two simulations may be required to get incoherent unpolarized result if the structure is not rotationally symmetric such as square array. That said, for square array, only one polarization is simulated.
• Angled incidence with broadband spectrum: The BFAST plane waves can be used with its own BCs for periodic structures. If it is for single wavelength, plane wave should be used with Bloch BCs.
• Boundary conditions: They are set in Boundary Conditions in FDTD Tab. Perpendicular to the direction of source injection axis PML BCs are always used to mimic light propagating infinitely. However in the transverse directions, boundary conditions must be matched with source type and structure properties. 

Monitors and analysis groups

For initial simulations a simple method to quantify the absorption can be used by only using transmission and reflection monitors. For heat simulation pab_adv analysis group can be used to generate the distributed absorption profile such as laser illumination. If the source is single wavelength, users can directly scale the real power when importing absorption into HEAT as the heat source;

The absorption analysis group from solar example has built_in function to output the absorption spatial density for one Sun. After multiplying the actual solar spectrum, the absorption spatial and spectral density has the unit in (W/m^3/nm). If it is broadband other than solar source, users will need to post process: having data for the real source spectrum and import it to the analysis group to replace solar spectrum, then obtain the real absorption density.

After simulation and analysis of the result such as absorption, refinement of the design is necessary. This can be done by parameter sweep and optimization. Parametrization is preferred for easy sweep and optimization, as well as for duplication. Such design may also be coupled with HEAT simulation, by separate simulations or by API.

Thermal simulations

Once optical simulation is completed, we can use thermal solver HEAT to simulate the temperature generated by the absorbed power. The steps are as follows:

Structure-devices

There are two ways to set up the structure: One is to set up objects from a completely new HEAT project by adding each individual geometries like what we have done in FDTD setup. The other way is to copy the objects or structure group from optical simulation by simply using Ctrl-C after selecting all the objects in FDTD project file, and then Ctrl-V to paste to HEAT.

Materials

The thermal material database is totally different from optical material database, therefore users will need to choose thermal materials from its material database. If users could not find the desired material in the database, user-created materials will need to be created.

Simulation region and mesh settings

Usually HEAT can use the FDTD's simulation sizes, except the substrate. Heat is conducted in solid so the substrate can be very thick, on the order of millimeters. To simplify the simulation, thermal impedance can be used in Boundary Conditions instead of the physical thickness.

The “min edge length” and the “max edge length” determines the grid properties, which can be swept for a converged results. In addition the NPs are very tiny so finer mesh may be needed. To better mesh curved surfaces the “deflection tolerance” in “Advanced Options” of “Mesh” can be modified from its default value.

If "mesh constraints" is used for a specific region to have finer mesh, it is the max edge length that is specified. Care must be taken not to use too small the max edge length.

Heat source

The heat source is mostly from optical simulation result and “import source” can be inserted in HEAT. Then import the result .mat file. The solar AM1.5 spectrum assumes 1000W/m^ of irradiance. For photocatalysis usually the power is much higher. The “scale factor” can be set for proper illumination.

Monitors

Usually “Temperature” monitor is enough to give result in the desired region. Additional "Power flow" monitor can be added if one wants to analyze the flow direction and power. In addition, from HEAT users can also visualize the result for the whole simulation region.

Boundary conditions

By default the HEAT solver provides insulation BCs. If other BCs such as fixed temperature are required, users will need to set and specify some parameters. In some cases thermal impedance can be specified to avoid large and uniform solid materials. The fluid material has built_in functionality that take care of its own boundary, no need to assign the solver BCs when fluid is at the boundary.

Material interfaces

HEAT solver solves thermal conduct in solid. So for solid materials their interfaces are automatically included. Air and other fluid materials have interfaces with solid to carry away some heat through convection or radiation. Thus the interface property must be specified at “Interface”. If the solid material is not in the list such as customer created-materials, users can manually add them to the interface with fluid materials. Then the thermal properties such as “Convection" and “Radiation” as well as required parameters can be set for each material that has interface with the fluid.

If material's thermal properties such as specific heat or thermal conductivity change significantly for the temperature range in simulation, users can turn on the thermal dependence with built_in models.

Once the set up finishes, temperature can be obtained after running the simulation. Not only one can get temperature and its distribution, the power flow at each boundary can also be extracted and analyzed by getting the data of "boundaries" in HEAT if necessary.

Since HEAT result depends on optical simulation, it is not necessary to change structure parameters however mesh settings may need to be tested to get converged results.

Based on the thermal simulation results, structure may be re-designed and redo the optical simulation in order to have the desired thermal results. For this purpose, co-simulation of HEAT together with FDTD can be performed using API.

Optical simulations

Users will need to design & optimize the plasmonic photocatalysis system to enhance the field intensity, increase the total absorption, reduce the reflection, and get the highest possible temperature etc by changing the NP geometry, arrangement and materials depending on the experiment goals.  This can be done through parameter sweep, optimization or API.  However, in this introductory example, we will only show how to set up optical simulation and get results for analysis and HEAT simulation.

Simulation settings

In this example, based on initial test, we choose the following parameters for illustration:

• NP Material: Au(Gold) - Palik
• Arrangement: Hexagonal array
• NP particles: sphere
• Sphere radius: 0.1um
• Period=0.21um
• Water layer thickness 0.1um to immerse the nano particles
• Substrate: SiO2 Palik

Above water, it is the background material of vacuum.

For simplicity, the water layer thickness is fixed to be the diameter of NPs. The NP array can be created manually, however for conciseness a structure group is set for NP radius, materials, period using script. The FDTD Xspan is the period, whereas the yspan is period*sqrt(3), and zspan is set to be 1um (about one wavelength). In x,y planes periodic BCs are used since the plane wave is at normal incidence. The rest are the default settings. With mesh accuracy 2, the small mesh size is about 15nm. Since the NPs are small in size, finer override mesh is added in the array region with mesh sizes of 5nm in x and 0.5*sqrt(3)/2 nm in y to have integer number of grids. Mesh size in z is also set to be 5nm.  Since the structure has symmetry, in "Advanced options" of "FDTD", select "force symmetric" in x/y.

Plane wave source with normal incidence towards -z in vacuum is used for this periodic structure. We do not extend to shorter wavelength than 0.4um of the solar spectrum in this example because it needs finer mesh. From initial test, longer wavelength from 0.8 um has small absorption. Thus wavelength range is set from 0.4 to 0.8 um. Since this array is hexagonal two perpendicular polarizations (E field along x with polarization angle 0, and along y with polarization angle 90 degree) are used in two separate simulations. To separate the absorption, output file name in the "absorption" analysis group will need to modify. However, for initial design, this group can be disabled and just optimize the absorption A=1-R-T where R and T are reflection and transmission. Further more, if the absorption of the two polarizations is about the same, again for initial test one polarization can be used to save time.

The built_in material fit has artificial resonance from Material Database in this wavelength range. For better description of Au material, the fitting parameters of "max coefficients" and "imaginary weight" in "show advanced" of the material explorer of Au (Gold) - Palik are set to be 7 and 10, As a result, the fitting RMS reduced to 0.13 from 1.3.

To obtain transmission and reflection, two frequency-domain power monitors are used.  After those settings, the set up is finished.

Results and discussion

One can run it using one polarization. If two polarizations are used, by setting AbsorptionEnable=0  in the script photocatalysis_2pol.lsf, the following results can be obtained. It can be seen that the simulated reflection, transmission and absorption of the two polarizations are very similar in this example.

To replicate the above results, download and run the script file, photocatalaysis_sim_2pol.lsf.

To get volumetric absorption density the "absorption" analysis group is used directly from Solar example. All the monitors record 50 frequency points set in"Monitors", "Global property".  The analysis group has redundant CHARGE results which can be removed by users. The useful output quantities from "absorption" analysis group are

• Pabs_total as a function of wavelength or frequency, which is the same as 1-R-T, normalized to source spectrum, dimensionless;
• Pabs: the raw absorption density as a function of space and wavelength, W/m^3. Its spatial integration is Pabs_total;
• Pabs_thermal: the spatial absorption density after integration over frequency, W/m^3.  
• Pabs_thermal_export is the same as Pabs_thermal except it contains specified multiple periods. In current simulation, only one period is required.

In the "absorption" analysis group, to count for semiconductor's bandgap, it has a parameter "bandgap", which is set to be 1.12ev, which corresponds a cut off frequency of 2.7e+014 Hz (1.1um). In this current example, no semiconductor material is involved so it is not necessary; in addition, the minimum frequency simulated is 3.33e+014 Hz, which is above the cutoff frequency. In other words, it does not affect the calculation of the absorption. Otherwise users will need to revise the script.  

To have the absorption data ready for HEAT use, users can enable "absorption" group, and assign "Q export filename" in a different name to output  Pabs_thermal_export in .mat format data used in HEAT. It saves time by use only one polarization for initial test.

For accurate analysis, two polarizations will be needed. In the script, setting AbsorptionEnable=1, and "Q_heat_x" and "Q_heat_y" are used. Once  the simulation is finished, two matlab data files are ready for HEAT simulation. The Figures below show the field intensity abs(E)^2 at different cross sections for incident source with y polarization.

It can be seen that the field is strongest at largest cross section of the NP array whereas at the surface it is smaller and at the bottom it is the smallest.

The figures below show the absorption density (W/m^3) at different cross sections (for visualization they are saturated in color) for incident source with y polarization.

The left figure shows the absorption is mainly in the upper portion of the spheres, and the strongest absorption is in the gaps of the spheres where the plasmonics plays roles. More detailed analysis can be carried for better understanding and designs.    

Result shows that for each polarization, the absorbed power in one period is about 0.013nW with direct Sun light, whereas the illumination power is about 0.076nW. The total absorption rate is around 34%.

Further design

The above structure is originated from a rough guess to illustrate the simulation setups and results. To have the optimal result, such as the highest absorption, users will need further work. For example, optimization and parameter sweep can target the absorption with one polarization since the other polarization gives similar result in this example, in order to save simulation and analysis time.

Tips

• Parameterization: The NP group is already parametrized including period, sphere radius, material. Other objects such as "water" and "substrate" can also be parametrized
• Due to symmetry of the geometry and the source polarization, symmetry BCs can be used to speed up the simulation since only one quarter simulation region is required.
• Unpolarized simulation: Source polarizations can be set in "model" for easy sweep and optimization, together with proper symmetry BCs, please refer this example: 3D pillar silicon solar cell. One can also set the proper symmetry BCs in the script file.
• Occasionally some users want to calculate the absorption cross section, which can be simply calculated as (1-R-T)*sourcepower/sourceinternsity.
• If the refractive index changes cannot be neglected due to high temperature users may need to modify it after thermal simulation and do iteration analysis.

Thermal simulations

Once optical simulation is completed, we can use  thermal solver HEAT to simulate the temperature generated by the absorbed power.

Simulation setup

Physical structures

We can copy the geometry from optical simulation and paste to HEAT, and then assign material properties from the material database:

Substrate material: “SiO2 (Glass) – Sze”

NP material: Set “mat” in the structure group to “Au (Gold) - CRC”. For better meshing, set the z position of this structure at 0.099um.

HEAT has no material data for water, so we need to find its thermal property and add a new material “water”. In nature water is fluid. In experiment, it is water that is heated up and vaporized at the boiling temperature and takes away the power dynamically. However in HEAT simulation of photocatalysis, the vaporization is not considered; and due to periodicity, no flow in the plane of array and thus no heat loss in theory. So we will set water material as an insulator with the following thermal parameters:

• density 996 kg/m^3
• specific heat 4178 J/kgK
• thermal conductivity 0.6W/mK

On top of water air is used as the fluid with its default thermal properties. In addition, the NP array is immersed inside water, so make sure the mesh order of "water" is 3, lower than the NPs.

In thermal simulation, the whole region must be fully filled with material, even it is air. For simplicity, we add a rectangle, name it as “Air”. Its dimensions are set to be large to cover the whole simulation region, but change its “mesh order” to 4. For better visualization, set is “alpha” in “Graphical rendering” to 0.2.

Simulation region and settings

Next, we add HEAT simulation region and choose 3D in "solver geometry". Set its x/y spans the same as the optical simulation, zmin is -2um (more discussions later), and zmax is 0.5um which is the same as optical simulation. Usually HEAT can use the FDTD's simulation sizes, except the substrate. Heat is conducted in solid so the substrate can be very thick, on the order of millimeters. To simplify the simulation, thermal impedance can be used in Boundary Conditions instead of the physical thickness.

The NPs are very tiny so finer mesh is needed. In addition, the simulation volume is small. So the “min edge length” is set 0.01um and “max edge length” is 0.1um in the HEAT “Mesh” settings.

Source and monitor

From the “Sources” section of the HEAT tab choose “Import heat”, an object called “heat” is created inside HEAT.  Set the import object “heat” at (0,0,0) and change its name to "heatX". Then import Q_heat_x.mat at “Data” and click “load”. In the “select attribute” choose “Pabs_thermal” to ensure the correct dataset. Then click “ok” to finish the import. Make sure the loaded “heat” object matches the geometry, without any spatial shift. If the other polarization is simulated, repeat the procedure to import Q_heat_y.mat and name it to be "heatY".

In this simulation, the absorption is 0.026nW which is much smaller than the typical power required for the photocatalysis(~ uW). So during the import process we set the “scale factor” to be 500 for each polarization (in total it is 1000 Sun), which could be realized by using a concentration lens.

A “Temperature” monitor is added for the array region in order to show temperature profile in this region.Usually “Temperature” monitor is enough to give result in the desired region. Additional "Power flow" monitor can be added if one wants to analyze the flow direction and power. In addition, from HEAT users can also visualize the result for the whole simulation region.

Boundary conditions and interfaces

The top (z max of HEAT) is air, which is a fluid material, the solver will take care its boundary condition.  The bottom of the substrate will be fixed to be room temperature by adding Boundary Conditions. Click to add "temperature" boundary in "Boundary Conditions" section of the HEAT tab, then edit it. Change its name to "bottom", and specifying it as “solver boundary” in z min; In the "General" tab set "bc mode" to be "steady state", set "sweep type" to be "single” at temperature of 300K.  If simulation substrate is thin, the heat created from nano particles will be quickly dissipated to room temperature since this boundary acts as a “cooler”. In addition, since we know the heat is mainly kept in the NP region, we can use thermal impedance to block the heat passing to the boundary while keeping the simulation length small.

The bulk thermal impedance is

$$R=\frac{L}{Ak}$$

where L is the length of the solid in meter, A is the cross-section in m², and k is the thermal conductivity in W/m-K.  Suppose a 1 mm thick substrate is used, for SiO2, the thermal impedance is about 9.5*10⁹ K/W where k=1.38 W/m-K. The thermal impedance inside the HEAT boundary condition is currently set to actually be the "Thermal insulance" which is R*A=0.72*10^-3 m²K/W.  Tick the "thermal impedance" and then input this value.

In addition, Air and any fluid material have interfaces with solid to carry away some heat through convection and radiation. We must specify the interface property. This is done at “boundary conditions” group by choosing the Air and the material with contact to Air. In this example such material is water, so add a convection boundary condition and choose "Air" and "water" as materials for the surface type "material:material". Also add a “Radiation” boundary condition between these two materials with its default parameter values. In addition, set the specific heat of Air in "Convection" to be h=10 W/m^2K.  Once the “Convection” is set the fluid will carry away some heat, and “Radiation” indicates that the fluid radiates way some heat, both are loss of power.  Set the emissivity of "water" in Air in the "Radiation" to be 0.9.

The simulation set up is complete.

Run and results

Then run the simulation. Once finished, from HEAT we can use the Visualizer and choose T for visualization. Local temperature can also be visualized using visualizer or exported using script from the "monitor":

aa=getresult("HEAT::monitor","temperature");
?T=max(aa.T);

A maximum temperature about 429K is obtained, which is much higher than the water boiling point of 373K. That said, the water has been already vaporized, which is not included in this simulation.

Local temperature distribution. Note, since the steady state has almost a uniform temperature in the NP region, the color is only for visualization purpose.  

The thermal simulation is complete.

Discussions

What is the possibly highest temperature for the given inject power of 0.26e^-12 W assuming it were a hot object with emissivity 0.9 ? When the environment temperature is 300K, this hot body would have a temperature about 1354K.  Next, test the highest possible temperature if insulation boundary condition at the bottom is used instead of fixed room temperature. By use of “insulating” in "model" of "bottom" BC, the only power loss is from the air since all other BCs are insulting. Simulation shows that the resulting temperature can be quite high, as high as about 1335K. This temperature would melt down the structure. In reality this cannot happen due to a number of reasons, such as the perfect insulation, water vaporization and  cold water flowing in etc. In such insulating case, the heat radiation plays important role and convection takes away only small amount of heat.

We suggest to use relatively coarse mesh for the initial simulations. Smaller min and max edge lengths may be used to get more accurate result, and "mesh constraint" can be used for regions which need finer meshing. However, due to possible discrepancy of the material thermal properties in simulation and experiment, mesh accuracy may not play important role for the final result.

Before you conclude the vaporization, a couple of things need to be sure:

• The thermal parameters can be quite different from their bulk values in nano scale; iteration simulation may be needed once compared with experiments;
• Insulating BCs are applied in the periodic directions however in experiment it can be difficult to have complete insulation;
• No cold water is injected into the simulation;
• When the temperature is high, the temperature dependence of Air thermal properties may need to be turned on, and convection parameters may need to match the experiment;
• If the temperature dependence of other materials' thermal properties are known, they can be included in the simulation.

Tips

• Since the structures are periodic in xy plane, in principle only one quarter can be simulated;
• to have a better geometry description of small curved surface such as sphere, “deflection tolerance” in “Advanced Options” of “Mesh” can be set to be smaller than the default 1nm;
• When large uniform materials are involved, the maximum edge length may need to increase;
• The imported source must be well resolved with finer meshing;  
• When the final temperature is predicted high, the "max update" temperature in HEAT-Advanced may be increased from its default 100K.

Related publications

1. Lin Zhou,Yingling Tan,Jingyang Wang,Weichao Xu,Ye Yuan,Wenshan Cai,Shining Zhu & Jia Zhu, 3D self-assembly of aluminium nanoparticles for plasmon-enhanced solar desalination,Nature Photonics 10, 393–398 (2016)
2. Oara Neumann, Alexander S. Urban, Jared Day, Surbhi Lal, Peter Nordlander, and Naomi J. Halas, Solar Vapor Generation Enabled by Nanoparticles, ACS Nano, 2013, 7 (1), pp 42–49.