Open tommyflynn13 opened 1 year ago
Hi Tommy,
That's a very interesting project, with lots of angles to investigate. I remember Goetzberger studied triangular geometries in one of his early papers, might be worth a look.
I don't know how amenable the LSC class is to extension to other geometries, but definitely worth looking in to. It was not something I considered when writing it.
The usual pvtrace trace approach is to, define a geometry, material and a light source and to throw rays into the screen and count where they come out.
I think you are on the right track starting with Mesh, later down the line you might want to write your own Triangle geometry class to handle the intersections yourself.
(Edit)
OK, so you have stuff more or less working. But probably I need to see the code to provide useful comments. If you have a fork, post a link to the line in question and the error is gives you.
No bother at all, I've been working through Jupyter and I'll attach the code below. Apologies if it's a bit messy, I'm picking it up on the fly!
Firstly, to create the polygon:
import trimesh
import shapely
from shapely.geometry import Point, Polygon
polygon = Polygon([[0, 0], [1, 2], [2, 0]]) # polygon
To extrude and visualise the polygon:
import logging
logging.getLogger('pvtrace').setLevel(logging.CRITICAL)
logging.getLogger('trimesh').setLevel(logging.CRITICAL)
from pvtrace import *
import trimesh
import shapely
# always need a world node
world = Node(
name="world",
geometry=Sphere(
radius=100.0
),
)
# Triangular Prism
mesh = Node(
name="mesh (triangular prism)",
geometry=Mesh(
trimesh=trimesh.creation.extrude_polygon(polygon, 0.5)
),
parent=world
)
mesh.translate((0.0, 0.0, 0.5))
scene = Scene(world)
vis = MeshcatRenderer(wireframe=True)
vis.render(scene)
vis.vis.jupyter_cell()
Finally, the adapted LSC.py code:
from pvtrace.material.component import Absorber, Luminophore
from pvtrace.light.light import Light, rectangular_mask
from pvtrace.light.event import Event
from pvtrace.scene.node import Node
from pvtrace.material.material import Material
from pvtrace.material.utils import isotropic, cone, lambertian
from pvtrace.scene.scene import Scene
from pvtrace.geometry.box import Box
from pvtrace.geometry.utils import EPS_ZERO
from pvtrace.data import lumogen_f_red_305, fluro_red
from pvtrace.scene.renderer import MeshcatRenderer
from pvtrace.material.surface import Surface, FresnelSurfaceDelegate
from pvtrace.material.distribution import Distribution
from pvtrace.algorithm import photon_tracer
from dataclasses import asdict
import numpy as np
import pandas as pd
import functools
import time
class OptionalMirrorAndSolarCell(FresnelSurfaceDelegate):
""" A delegate adds an ideal specular mirror to the bottom surface and
perfectly indexed matched and perfectly absorbing solar cells to the edges.
"""
def __init__(self, lsc):
super(OptionalMirrorAndSolarCell, self).__init__()
self.lsc = lsc
def reflectivity(self, surface, ray, geometry, container, adjacent):
normal = geometry.normal(ray.position)
cell_locations = self.lsc._solar_cell_surfaces
back_surface_mirror = self.lsc._back_surface_mirror_info
want_back_surface_mirror = back_surface_mirror["want_back_surface_mirror"]
if np.allclose((0, 0, -1), normal) and want_back_surface_mirror:
return 1.0 # perfect mirror
elif np.allclose((-1, 0, 0), normal) and "left" in cell_locations: # left
return 0.0 # perfect absorption
elif np.allclose((1, 0, 0), normal) and "right" in cell_locations: # right
return 0.0
elif np.allclose((0, -1, 0), normal) and "near" in cell_locations: # near
return 0.0
elif np.allclose((0, 1, 0), normal) and "far" in cell_locations: # far
return 0.0
return super(OptionalMirrorAndSolarCell, self).reflectivity(
surface, ray, geometry, container, adjacent
) # opt-out of handling custom reflection
def transmitted_direction(self, surface, ray, geometry, container, adjacent):
cell_locations = self.lsc._solar_cell_surfaces
normal = geometry.normal(ray.position)
if (
(np.allclose((-1, 0, 0), normal) and "left" in cell_locations)
or (np.allclose((1, 0, 0), normal) and "right" in cell_locations)
or (np.allclose((0, -1, 0), normal) and "near" in cell_locations)
or (np.allclose((0, 1, 0), normal) and "far" in cell_locations)
):
return ray.direction # solar cell is perfectly index matched
return super(OptionalMirrorAndSolarCell, self).transmitted_direction(
surface, ray, geometry, container, adjacent
) # opt-out of handling custom reflection
class AirGapMirror(FresnelSurfaceDelegate):
def __init__(self, lsc):
super(AirGapMirror, self).__init__()
self.lsc = lsc
def reflectivity(self, surface, ray, geometry, container, adjacent):
return 1.0 # perfect reflector. NB don't reduce.
def reflected_direction(self, surface, ray, geometry, container, adjacent):
if not self.lsc._air_gap_mirror_info["lambertian"]:
# Specular reflection
return super(AirGapMirror, self).transmitted_direction(
surface, ray, geometry, container, adjacent
)
else:
normal = geometry.normal(ray.position)
if not np.allclose((0.0, 0.0, 1.0), normal): # top surface
raise NotImplementedError("Not yet generalised to other surfaces.")
# Currently this return lambertian direction along +z axis and is not
# generalised to other orientations. This is simple to do using a transform
# which first moves into to the z+ frame and then back out.
return tuple(lambertian().tolist())
class LSC(object):
"""Abstraction of a luminescent solar concentrator.
This is intended to be a high-level API to easy use.
"""
def __init__(self, wavelength_range=None, n0=1.0, n1=1.5):
super(LSC, self).__init__()
if wavelength_range is None:
self.wavelength_range = np.arange(400, 800)
self.n0 = n0
self.n1 = n1
self._solar_cell_surfaces = set()
self._back_surface_mirror_info = {"want_back_surface_mirror": False}
self._air_gap_mirror_info = {"want_air_gap_mirror": False, "lambertian": False}
self._scene = None
self._renderer = None
self._store = None
self._df = None
self._counts = None
self._user_lights = []
self._user_components = []
def _make_default_components(self):
""" Default LSC contains Lumogen F Red 305. With concentration such that
the absorption coefficient at peak is 10 cm-1.
"""
x = self.wavelength_range
coefficient = lumogen_f_red_305.absorption(x) * 10.0 # cm-1
emission = lumogen_f_red_305.emission(x)
coefficient = np.column_stack((x, coefficient))
emission = np.column_stack((x, emission))
lumogen = {
"cls": Luminophore,
"name": "Lumogen F Red 305",
"coefficient": coefficient,
"emission": emission,
"quantum_yield": 1.0,
"phase_function": None, # will select isotropic
}
background = {"cls": Absorber, "coefficient": 0.1, "name": "Background"} # cm-1
return [lumogen, background]
def _make_default_lights(self):
""" Default light is a spotlight (cone of 20-deg) of single wavelength 555nm.
"""
light = {
"name": "Light",
"location": (0.0, 0.0, self.size[-1] * 5),
"rotation": (np.radians(180), (1, 0, 0)),
"direction": functools.partial(cone, np.radians(20)),
"wavelength": None,
"position": None,
}
return [light]
def _make_scene(self):
""" Creates the scene based on configuration values.
"""
# Make world node
world = Node(
name="World",
geometry=Box(
(100, 100, 100), material=Material(refractive_index=self.n0)
),
)
# Create components (Absorbers, Luminophores and Scatteres)
if len(self._user_components) == 0:
self._user_components = self._make_default_components()
components = []
for component_data in self._user_components:
cls = component_data.pop("cls")
coefficient = component_data.pop("coefficient")
component = cls(coefficient, **component_data)
components.append(component)
# Create LSC node
lsc = Node(
name="LSC",
geometry=Mesh(
trimesh=trimesh.creation.extrude_polygon(polygon, 0.5),
material=Material(
refractive_index=self.n1,
components=components,
surface=Surface(delegate=OptionalMirrorAndSolarCell(self)),
),
),
parent=world,
)
if self._air_gap_mirror_info["want_air_gap_mirror"]:
sheet_thickness = 0.1 # make it appear thinner than the LSC
air_gap_mirror = Node(
name="Air Gap Mirror",
geometry=Mesh(
trimesh=trimesh.creation.extrude_polygon(polygon, 0.1), # same surface air but very thin
material=Material(
refractive_index=self.n0,
components=[],
surface=Surface(delegate=AirGapMirror(self)),
),
),
parent=world,
)
# Move adjacent to bottom surface with a small air gap
air_gap_mirror.translate((0.0, 0.0, -(0.5 * 0.5 + 0.1)))
# Use user light if any have been given, otherwise use default values.
if len(self._user_lights) == 0:
self._user_lights = self._make_default_lights()
# Create light nodes
for light_data in self._user_lights:
name = light_data["name"]
light = Light(
name=name,
direction=light_data["direction"],
wavelength=light_data["wavelength"],
position=light_data["position"],
)
light_node = Node(name=name, light=light, parent=world)
light_node.location = light_data["location"]
if light_data["rotation"]:
light_node.rotate(*light_data["rotation"])
self._scene = Scene(world)
def component_names(self):
if self._scene is None:
raise ValueError("Run a simulation before calling this method.")
return {c["name"] for c in self._user_components}
def light_names(self):
if self._scene is None:
raise ValueError("Run a simulation before calling this method.")
return {l["name"] for l in self._user_lights}
def add_luminophore(
self, name, coefficient, emission, quantum_yield, phase_function=None
):
self._user_components.append(
{
"cls": Luminophore,
"name": name,
"coefficient": coefficient,
"emission": emission,
"quantum_yield": quantum_yield,
"phase_function": phase_function,
}
)
def add_absorber(self, name, coefficient):
self._user_components.append(
{"cls": Absorber, "name": name, "coefficient": coefficient}
)
def add_scatterer(self, name, coefficient, phase_function=None):
self._user_components.append(
{
"cls": Scatterer,
"name": name,
"coefficient": coefficient,
"phase_function": phase_function,
}
)
def add_light(
self,
name,
location, # node location in parent
rotation=None, # node rotation in parent frame
direction=None, # direction delegate callable
wavelength=None, # wavelength delegate callable
position=None, # position delegate callable
):
self._user_lights.append(
{
"name": name,
"location": location,
"rotation": rotation,
"direction": direction,
"wavelength": wavelength,
"position": position,
}
)
def add_solar_cell(self, facets):
if not isinstance(facets, (list, tuple, set)):
raise ValueError("Facets should be a set. e.g. `{'left', 'right'}`")
facets = set(facets)
allowed = {"left", "near", "far", "right"}
if not facets.issubset(allowed):
raise ValueError("Solar cell have allowed surfaces", allowed)
self._solar_cell_surfaces = facets.union(self._solar_cell_surfaces)
def add_back_surface_mirror(self):
self._back_surface_mirror_info = {"want_back_surface_mirror": True}
def add_air_gap_mirror(self, lambertian=False):
self._air_gap_mirror_info = {
"want_air_gap_mirror": True,
"lambertian": lambertian,
}
# Simulate
def show(
self,
wireframe=True,
baubles=True,
bauble_radius=None,
world_segment="short",
short_length=None,
open_browser=False,
):
if bauble_radius is None:
bauble_radius = np.min(self.size) * 0.05
if short_length is None:
short_length = np.min(self.size) * 0.1
self._add_history_kwargs = {
"bauble_radius": bauble_radius,
"baubles": baubles,
"world_segment": world_segment,
"short_length": short_length,
}
if self._scene is None:
self._make_scene()
self._renderer = MeshcatRenderer(
open_browser=open_browser,
transparency=False,
opacity=0.5,
wireframe=wireframe,
max_histories=50,
)
self._renderer.render(self._scene)
time.sleep(1.0)
return self._renderer
def simulate(self, n, progress=None, emit_method="kT"):
if self._scene is None:
self._make_scene()
scene = self._scene
# Simulate can be called multiple time to append rays to the store
if self._store is None:
store = {"entrance_rays": [], "exit_rays": []}
vis = self._renderer
count = 0
for ray in scene.emit(n):
history = photon_tracer.follow(scene, ray, emit_method=emit_method)
rays, events = zip(*history)
store["entrance_rays"].append((rays[1], events[1]))
if events[-1] in (Event.ABSORB, Event.KILL):
# final event is a lost store path information at final event
store["exit_rays"].append((rays[-1], events[-1]))
elif events[-1] == Event.EXIT:
# final event hits the world node. Store path information at
# penultimate location
store["exit_rays"].append((rays[-2], events[-2]))
# Update visualiser
if vis:
vis.add_history(history, **self._add_history_kwargs)
# Progress callback
if progress:
count += 1
progress(count)
self._store = store
print("Tracing finished.")
print("Preparing results.")
df = self._make_dataframe()
df = self.expand_coords(df, "direction")
df = self.expand_coords(df, "position")
df = self.label_facets(df, *self.size)
self._df = df
def _make_dataframe(self):
# to-do: Only need to process additional rays not whole dataframe! Optimise!
df = pd.DataFrame()
# Rays entering the scene
for ray, event in self._store["entrance_rays"]:
rep = asdict(ray)
rep["kind"] = "entrance"
rep["event"] = event.name.lower()
df = df.append(rep, ignore_index=True)
# Rays exiting the scene
for ray, event in self._store["exit_rays"]:
rep = asdict(ray)
rep["kind"] = "exit"
rep["event"] = event.name.lower()
df = df.append(rep, ignore_index=True)
self._df = df
return df
def expand_coords(self, df, column):
""" Returns a dataframe with coordinate column expanded into components.
Parameters
----------
df : pandas.DataFrame
The dataframe
column : str
The column label
Returns
-------
df : pandas.DataFrame
The dataframe with the column expanded.
Example
-------
Given the dataframe::
df = pd.DataFrame({'position': [(1,2,3)]})
the function will return a new dataframe::
edf = expand_coords(df, 'position')
edf == pd.DataFrame({'position_x': [1], 'position_y': [2], 'position_z': [3]})
"""
coords = np.stack(df[column].values)
df["{}_x".format(column)] = coords[:, 0]
df["{}_y".format(column)] = coords[:, 1]
df["{}_z".format(column)] = coords[:, 2]
df.drop(columns=column, inplace=True)
return df
def label_facets(self, df, length, width, height):
""" Label rows with facet names for a box LSC.
Notes
-----
This function only works if the coordinates in the dataframe
are in the local frame of the box. If the coordinates are in the
world frame then this will still work provided the box is axis
aligned with the world node and centred at the origin.
"""
xmin, xmax = -0.5 * length, 0.5 * length
ymin, ymax = -0.5 * width, 0.5 * width
zmin, zmax = -0.5 * height, 0.5 * height
df.loc[(np.isclose(df["position_x"], xmin, atol=EPS_ZERO)), "facet"] = "left"
df.loc[(np.isclose(df["position_x"], xmax, atol=EPS_ZERO)), "facet"] = "right"
df.loc[(np.isclose(df["position_y"], ymin, atol=EPS_ZERO)), "facet"] = "far"
df.loc[(np.isclose(df["position_y"], ymax, atol=EPS_ZERO)), "facet"] = "near"
df.loc[(np.isclose(df["position_z"], zmin, atol=EPS_ZERO)), "facet"] = "bottom"
df.loc[(np.isclose(df["position_z"], zmax, atol=EPS_ZERO)), "facet"] = "top"
return df
def _make_counts(self, df):
if self._counts is not None:
return self._counts
components = self._scene.component_nodes
lights = self._scene.light_nodes
all_components = {component.name for component in components}
all_lights = {light.name for light in lights}
# Count solar photons exiting
solar_out = dict()
for facet in {"left", "right", "near", "far", "top", "bottom"}:
solar_out[facet] = self.spectrum(
facets={facet}, source=all_lights, kind="last"
).shape[0]
# Count solar photons entering
solar_in = dict()
for facet in {"left", "right", "near", "far", "top", "bottom"}:
solar_in[facet] = self.spectrum(
facets={facet}, source=all_lights, kind="first"
).shape[0]
# Count luminescent photons exiting
lum_out = dict()
for facet in {"left", "right", "near", "far", "top", "bottom"}:
lum_out[facet] = self.spectrum(
facets={facet}, source=all_components, kind="last"
).shape[0]
# Count luminescent photons entering
lum_in = dict()
for facet in {"left", "right", "near", "far", "top", "bottom"}:
lum_in[facet] = self.spectrum(
facets={facet}, source=all_components, kind="first"
).shape[0]
self._counts = counts = pd.DataFrame(
{
"Solar In": pd.Series(solar_in),
"Solar Out": pd.Series(solar_out),
"Luminescent Out": pd.Series(lum_out),
"Luminescent In": pd.Series(lum_in),
},
index=["left", "right", "near", "far", "top", "bottom"],
)
return counts
def spectrum(self, facets=set(), kind="last", source="all", events=None):
if self._df is None:
raise ValueError("Run a simulation before calling this method.")
df = self._df
if kind is not None:
if not kind in {"first", "last"}:
raise ValueError("Direction must be either `'first'` or `'last'.`")
if kind is None:
want_kind = True # Opt-out
else:
if kind == "first":
want_kind = df["kind"] == "entrance"
else:
want_kind = df["kind"] == "exit"
all_sources = self.component_names() | self.light_names()
if source == "all":
want_sources = all_sources
else:
if isinstance(source, str):
source = {source}
if not set(source).issubset(all_sources):
unknown_source_set = set(source).difference(all_sources)
raise ValueError("Unknown source requested.", unknown_source_set)
if source == "all":
want_source = df["source"].isin(all_sources)
else:
want_source = df["source"].isin(set(source))
if isinstance(facets, (list, tuple, set)):
if len(facets) > 0:
want_facets = df["facet"].isin(set(facets))
else:
want_facets = True # don't filter by facets
else:
raise ValueError(
"`'facets'` should be a set `{'left', 'right'}`", {"got": facets}
)
if events is None:
want_events = True # Don't filter by events
else:
all_events = {e.name.lower() for e in Event}
if isinstance(events, (list, tuple, set)):
events = set(events)
if events.issubset(all_events):
want_events = df["event"].isin(events)
else:
raise ValueError(
"Contained some unknown events",
{"got": events, "expected": all_events},
)
else:
raise ValueError(
"Events must be set of event strings", {"allowed": all_events}
)
return df.loc[want_kind & want_source & want_facets & want_events]["wavelength"]
def counts(self):
df = self._df
if df is None:
df = self._make_dataframe()
df = self.expand_coords(df, "direction")
df = self.expand_coords(df, "position")
df = self.label_facets(df, *self.size)
counts = self._make_counts(df)
return counts
def summary(self):
counts = self._make_counts(self._df)
all_facets = {"left", "right", "near", "far", "top", "bottom"}
lum_collected = 0
for facet in self._solar_cell_surfaces:
lum_collected += counts["Luminescent Out"][facet]
lum_escaped = 0
for facet in all_facets - self._solar_cell_surfaces:
lum_escaped += counts["Luminescent Out"][facet]
incident = 0
for facet in all_facets:
incident += counts["Solar In"][facet]
lost = self.spectrum(source="all", events={"absorb"}, kind="last").shape[0]
optical_efficiency = lum_collected / incident
waveguide_efficiency = lum_collected / (lum_collected + lum_escaped)
(l, w, d) = self.size
a1 = w * l
a2 = 2 * l * d + 2 * w * d
Cg = a1 / a2
n = self.n1
s = pd.Series(
{
"Optical Efficiency": optical_efficiency,
"Waveguide Efficiency": waveguide_efficiency,
"Waveguide Efficiency (Thermodynamic Prediction)": (
n ** 2 / (Cg + n ** 2)
),
"Non-radiative Loss (fraction):": lost / incident,
"Incident": incident,
"Geometric Concentration": Cg,
"Refractive Index": n,
"Cell Surfaces": self._solar_cell_surfaces,
"Components": self.component_names(),
"Lights": self.light_names(),
}
)
return s
def report(self):
print()
print("Simulation Report")
print("-----------------")
print()
print("Surface Counts:")
print(self.counts())
print()
print("Summary:")
print(self.summary())
The main changes made to the LSC.py script were in the LSC() class. I removed the size parameter in definit, updated the world node as it involved a multiple of the (l, w, d) parameters which no longer exist, and updated the geometry of the LSC node to trimesh.creation.extrude_polygon instead of Box().
lsc = LSC((polygon, 1.0)) # size in cm
lsc.show() # open visualiser
lsc._renderer.vis.jupyter_cell()
The errors appear when I run the above script to run the LSC report. The errors mention things like the self.size parameter in the bauble_radius code, which I assume means it is attached to the size of the LSC box which no longer exists. I've tried tinkering with this but my python knowledge is very limited so I'm not sure where to begin.
Hi Tommy,
I think this an error from the MeshcatRenderer. Baubles are the little spheres that appear on the ray paths when it hits something, they look like Christmas tree decorations hence the (bad) name.
Can you put the actual traceback here?
But I think the problem is that you are, perhaps, giving pvtrace the trimesh mesh object. You will need to wrap the actual trimesh polygon object with using pvtrace's Mesh class (https://github.com/danieljfarrell/pvtrace/blob/master/pvtrace/geometry/mesh.py). This implement the Geometry protocol so that pvtrace and query the geometry object.
So rather than
LSC((polygon, 1.0))
try,
LSC((Mesh(polygon), 1.0))
No problem at all, I'll past the traceback here below:
AttributeError Traceback (most recent call last)
Input In [7], in <cell line: 2>()
1 lsc = LSC((polygon, 1.0)) # size in cm
----> 2 lsc.show() # open visualiser
3 lsc._renderer.vis.jupyter_cell()
Input In [6], in LSC.show(self, wireframe, baubles, bauble_radius, world_segment, short_length, open_browser)
299 def show(
300 self,
301 wireframe=True,
(...)
306 open_browser=False,
307 ):
309 if bauble_radius is None:
--> 310 bauble_radius = np.min(self.size) * 0.05
312 if short_length is None:
313 short_length = np.min(self.size) * 0.1
AttributeError: 'LSC' object has no attribute 'size'
And if you change line 1 to be
LSC((Mesh(polygon), 1.0))
Then it gives this:
---------------------------------------------------------------------------
AttributeError Traceback (most recent call last)
Input In [8], in <cell line: 1>()
----> 1 lsc = LSC((Mesh(polygon), 1.0)) # size in cm
2 lsc.show() # open visualiser
3 lsc._renderer.vis.jupyter_cell()
File ~\miniconda3\envs\pvtrace-env\lib\site-packages\pvtrace\geometry\mesh.py:17, in Mesh.__init__(self, trimesh, material)
15 def __init__(self, trimesh, material=None):
16 super(Mesh, self).__init__()
---> 17 trimesh.vertices -= trimesh.center_mass
18 self.trimesh = trimesh
19 self._material = material
AttributeError: 'Polygon' object has no attribute 'vertices'
Edit: I also tried this input:
lsc = LSC(Mesh(trimesh=trimesh.creation.extrude_polygon(polygon, 1.0))) # size in cm
lsc.show() # open visualiser
lsc._renderer.vis.jupyter_cell()
Which gave the traceback:
---------------------------------------------------------------------------
AttributeError Traceback (most recent call last)
Input In [11], in <cell line: 2>()
1 lsc = LSC(Mesh(trimesh=trimesh.creation.extrude_polygon(polygon, 1.0))) # size in cm
----> 2 lsc.show() # open visualiser
3 lsc._renderer.vis.jupyter_cell()
Input In [6], in LSC.show(self, wireframe, baubles, bauble_radius, world_segment, short_length, open_browser)
299 def show(
300 self,
301 wireframe=True,
(...)
306 open_browser=False,
307 ):
309 if bauble_radius is None:
--> 310 bauble_radius = np.min(self.size) * 0.05
312 if short_length is None:
313 short_length = np.min(self.size) * 0.1
AttributeError: 'LSC' object has no attribute 'size'
What's the type of the polygon
?
type(polygon)
It should be a trimesh.Trimesh
or subclass.
It says shapely.geometry.polygon.Polygon
🤔 try and convert that to a trimesh.Trimesh object and it "should just work".
No bother, I'll give that a go and let you know how I get on! Should I keep using the LSC.py I updated and removed the 'size' parameters, or would I be best to revert back to your original LSC.py script?
Hi Daniel,
Just to give you an update I tried a new approach and it seems to be coming together! Updated the LSC node script to be
# Create LSC node
polygon = Polygon([[0, 0], [(w/2), l], [w, 0]])
lsc = Node(
name="LSC",
geometry=Mesh(
trimesh=trimesh.creation.extrude_polygon(polygon, d),
material=Material(
refractive_index=self.n1,
components=components,
surface=Surface(delegate=OptionalMirrorAndSolarCell(self)),
),
),
parent=world,
)
I've pasted what appeared in the visualiser after running the new LSC.py script
.
Oh looks like it's doing something. Think you solved the initial problem. Nice!
After updating the air gap to be triangular as well, the visualiser worked however the report gave a new traceback error:
GeometryError Traceback (most recent call last)
Input In [20], in <cell line: 754>()
751 lsc.report()
754 if __name__ == "__main__":
--> 755 test2()
Input In [20], in test2()
748 lsc.show()
749 throw = 300
--> 750 lsc.simulate(throw, emit_method="redshift")
751 lsc.report()
Input In [20], in LSC.simulate(self, n, progress, emit_method)
350 count = 0
351 for ray in scene.emit(n):
--> 352 history = photon_tracer.follow(scene, ray, emit_method=emit_method)
353 rays, events = zip(*history)
354 store["entrance_rays"].append((rays[1], events[1]))
File ~\miniconda3\envs\pvtrace-env\lib\site-packages\pvtrace\algorithm\photon_tracer.py:194, in follow(scene, ray, maxsteps, maxpathlength, emit_method)
192 ray = ray.representation(scene.root, hit)
193 if surface.is_reflected(ray, hit.geometry, container, adjacent):
--> 194 ray = surface.reflect(ray, hit.geometry, container, adjacent)
195 ray = ray.representation(hit, scene.root)
196 history.append((ray, Event.REFLECT))
File ~\miniconda3\envs\pvtrace-env\lib\site-packages\pvtrace\material\surface.py:245, in Surface.reflect(self, ray, geometry, container, adjacent)
242 def reflect(self, ray, geometry, container, adjacent):
243 """ Returns ray which is reflected from the interface.
244 """
--> 245 direction = self.delegate.reflected_direction(
246 self, ray, geometry, container, adjacent
247 )
248 if not isinstance(direction, tuple):
249 raise ValueError(
250 "Delegate method `reflected_direction` should return a tuple."
251 )
Input In [20], in AirGapMirror.reflected_direction(self, surface, ray, geometry, container, adjacent)
74 def reflected_direction(self, surface, ray, geometry, container, adjacent):
75 if not self.lsc._air_gap_mirror_info["lambertian"]:
76 # Specular reflection
---> 77 return super(AirGapMirror, self).transmitted_direction(
78 surface, ray, geometry, container, adjacent
79 )
80 else:
81 normal = geometry.normal(ray.position)
File ~\miniconda3\envs\pvtrace-env\lib\site-packages\pvtrace\material\surface.py:173, in FresnelSurfaceDelegate.transmitted_direction(self, surface, ray, geometry, container, adjacent)
171 n2 = adjacent.geometry.material.refractive_index
172 # Be tolerance with definition of surface normal
--> 173 normal = geometry.normal(ray.position)
174 if np.dot(normal, ray.direction) < 0.0:
175 normal = flip(normal)
File ~\miniconda3\envs\pvtrace-env\lib\site-packages\pvtrace\geometry\mesh.py:75, in Mesh.normal(self, surface_point)
71 raise GeometryError(
72 "Mesh must have a single closest point to calculate normal."
73 )
74 if not np.any(np.absolute(distances) < EPS_ZERO):
---> 75 raise GeometryError(
76 "Point is not on surface.",
77 {
78 "point": surface_point,
79 "geometry": self,
80 "distances": distances,
81 "threshold": EPS_ZERO,
82 },
83 )
84 normal = tuple(mesh.face_normals[triangle_id[0]])
85 return normal
GeometryError: ('Point is not on surface.', {'point': (-0.48078641173750625, 0.635805804066633, -0.006357771365071763), 'geometry': <pvtrace.geometry.mesh.Mesh object at 0x000001581DF88BE0>, 'distances': array([6.07280425e-13]), 'threshold': 2.220446049250313e-13})
This seems to be a common issue with Mesh geometry. It's related to the precision of the intersection points.
Cab you define your scene in cm or mm rather than meters? This will make everything bigger so thresholding needs less precision.
Everything seems to be in working order now, I've attached an image sjowing a quick simulation with 100 rays. My (hopefully!) last question, I was just wondering would you have any idea how I would alter the surfaces in the simulation report? Currently they seem to be (Left, right, top, bottom) as it would be for a box geometry, I imagine I should change the coordinates to have just (left, right, bottom) for mine, with the left and right both being sloped surfaces rather than straight as they would have been in the box.
Simulation Report
-----------------
Surface Counts:
Solar In Solar Out Luminescent Out Luminescent In
left 0 0 0 0
right 0 0 0 0
near 0 0 0 0
far 0 0 0 0
top 100 3 16 0
bottom 0 1 14 0
Summary:
Optical Efficiency 0.0
Waveguide Efficiency 0.0
Waveguide Efficiency (Thermodynamic Prediction) 0.473684
Non-radiative Loss (fraction): 0.28
Incident 100
Geometric Concentration 2.5
Refractive Index 1.5
Cell Surfaces {}
Components {Lumogen F Red 305, Background}
Lights {Light}
dtype: object
Thanks a million for all your help this far also!
Take a look at label_facets
function in LSC.py
.
This assigns labels to the points sitting on the surface. For example, if x-coordinate of the of the point is on the box plane with x=xmin we definite that as being the left surface.
Remember you are comparing floats so need to have a tolerance.
You will have to come up with a scheme for doing this labelling with your geometry.
It will be the same for the edge surfaces and top and bottom, but they hypotenuse will need more of a complicated comparison to evaluate.
I've taken a look at that piece of code, and updated it to the following. Not sure if it makes sense as I don't recognise the formula or attributes. just been playing around with it until it ran successfully.
def label_facets(self, df, length, width, height):
""" Label rows with facet names for a box LSC.
Notes
-----
This function only works if the coordinates in the dataframe
are in the local frame of the box. If the coordinates are in the
world frame then this will still work provided the box is axis
aligned with the world node and centred at the origin.
"""
xmin, xmax = -0.5 * width, 0.5 * width
ymin, ymax = -0.3334 * length, 0.6666 * length
zmin, zmax = -0.5 * height, 0.5 * height
df.loc[((df["position_x"] < (width/2))), "facet"] = "left"
df.loc[((df["position_x"] > (width/2))), "facet"] = "right"
df.loc[(np.isclose(df["position_y"], ymin, atol=EPS_ZERO)), "facet"] = "near"
df.loc[(np.isclose(df["position_z"], zmin, atol=EPS_ZERO)), "facet"] = "bottom"
df.loc[(np.isclose(df["position_z"], zmax, atol=EPS_ZERO)), "facet"] = "top"
return df
The output was as follows:
Simulation Report
-----------------
Surface Counts:
Solar In Solar Out Luminescent Out Luminescent In
left 127 163 34 0
right 89 89 0 0
near 0 0 0 0
top 84 2 12 0
bottom 0 0 0 0
Seems to be giving "Solar In" values for the left and right surfaces, so something must be incorrect in the way I currently have it written for the "left" and "right" facet labels.
I'd suggest continuing with your notebook approach, as this it probably a bit too complicated for now.
Look at the the Quick Start notebook and others in the examples folder. They show how to throw a ray into the scene and record where it escapes.
It too early to make these neat for your geometry. Just get it to work first and then worry about making it simple to use.
The LSC.py module was actually the last thing I coded!
Assigning positions to facets it actually quite tricky if just relying on the coordinates. I think pvtrace should probably make this easier by retuning an ID.
There is a refactor in the works to improve this, but I've been a bit busy the last few years 🥴
You might be right! Just giving it a final shot as I feel I could be close enough. At the minute, I'm attmepting to define the equation of the line (y1 and y2) on the left and right hand sides, and defining "left" and "right" as being close to that line.
def label_facets(self, df, length, width, height):
""" Label rows with facet names for a box LSC.
Notes
-----
This function only works if the coordinates in the dataframe
are in the local frame of the box. If the coordinates are in the
world frame then this will still work provided the box is axis
aligned with the world node and centred at the origin.
"""
y1 = ((width/length) * df["position_x"]) + ((2/3) * length)
y2 = -((width/length) * df["position_x"]) + ((2/3) * length)
xmin, xmax = -0.5 * width, 0.5 * width
ymin, ymax = -(1/3) * length, (2/3) * length
zmin, zmax = -0.5 * height, 0.5 * height
df.loc[(np.isclose(df["position_y"], y1, atol=EPS_ZERO)), "facet"] = "left"
df.loc[(np.isclose(df["position_y"], y2, atol=EPS_ZERO)), "facet"] = "right"
df.loc[(np.isclose(df["position_y"], ymin, atol=EPS_ZERO)), "facet"] = "near"
df.loc[(np.isclose(df["position_y"], ymax, atol=EPS_ZERO)), "facet"] = "far"
df.loc[(np.isclose(df["position_z"], zmin, atol=EPS_ZERO)), "facet"] = "bottom"
df.loc[(np.isclose(df["position_z"], zmax, atol=EPS_ZERO)), "facet"] = "top"
return df
It doesn't seem to be working as the results aren't adding up, but do you think I could be on the right track with this? Would I be better off trying to define an array of points along the line (the equation of the line is y = 1.4x + 4.666 for the dimensions of my LSC) and using that, and if so how would I do that?
I've attached an image from Desmos of the equations of the three lines which make up my triangular geometry.
Thanks again, Tommy.
Hi Daniel,
Firstly just wanted to say thanks a million for the work you've done on this code, the software is very impressive and easy to follow for a beginner like myself! I was shown your code by a fellow university student, and am hoping to manipulate it to simulate a triangular LSC panel being fabricated for use in a geodesic dome installation as part of my Masters dissertation.
Not an issue as such, but I was hoping to pick your brain as to how you would edit the LSC.py code to work with a triangular prism rather than a box? As I mentioned I am a complete beginner with Python, I have figured out how to extrude a triangular polygon using shapely and trimesh, and have passed light through it using one of your tutorials. I am struggling however to fully adapt the LSC.py code to work with this mesh geometry, specifically as the "size" and "self.size" parameters don't seem to apply to the (polygon, height) input using trimesh.creation.extrude_polygon.
Any help would be greatly appreciated.
Regards, Thomas.